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NAVY  INTERIM  SURFACE  SHIP  MODEL  (NISSM)  II 


INTRODUCTION 

The  Navy  Interim,  Surface  Ship  Model  (NISSM)  II  is  a  FORTRAN  V  computer 
program  designed  to  predict  the  performance  of  existing  and  proposed  active 
sonar  systems  using  ray-tracing  techniques.  The  measure  of  performance  is 
expressed  in  terms  of  probability  of  detection  versus  target  range  for  a  given 
false-alarm  rate.  Several  intermediate  computations,  including  ray  trajectories, 
propagation  losses,  and  boundary  and  volume  reverberation,  are  also  displayed. 

Among  the  noteworthy  features  of  NISSM  II  is  the  manner  in  which  the  velocity 
of  sound  in  the  ocean  is  represented.  Whereas  most  ray-tracing  techniques  use 
a  piecewise  constant  gradient  fit,  NISSM  II  approximates  the  sound  velocity  with 
a  continuous  function  of  depth  having  a  continuous  gradient.  This  reduces  the 
problem  of  false  caustics,  a  computational  difficulty  that  may  be  encountered  in 
the  prediction  of  geometrical  spreading  loss.  Since  classical  ray  tracing  is  a 
high-frequency  solution  of  the  wave  equation  and  does  not  accurately  describe 
low-frequency  diffraction  phenomena,  the  ray  theory  is  modified  to  improve 
results  (pp.  21-26). 

The  most  time  consuming  phase  in  executing  NISSM  II  is  the  volume- 
reverberation  computation.  The  volume-scattering  strength,  a  function  of  depth, 
is  expressed  in  terms  of  scattering  strength  per  unit  volume— not  as  scattering 
strength  per  unit  area  of  water  column.  A  careful  integration  of  the  volume- 
reverberation  integral  is  performed,  whereby  the  complicated  multipath  structures 
that  maybe  encountered  in  practice  are  considered.  Typical  volume- reverberation 
computations  by  the  Univac  1108  Multiprocessor  System  require  3  min.  A  plot 
of  probability  of  detection  versus  target  range  requires  approximately  3. 5  min. 

On  the  other  hand,  NISSM  II  can  produce  a  ray  diagram  or  a  propagation-loss 
curve  in  10  to  20  s. 

The  major  limitations  of  the  NISSM  II  model  are  that  (1)  the  ocean  bottom 
must  be  horizontal,  (2)  the  velocity  of  sound  must  not  vary  with  range,  and  (3)  the 
usual  constraints  on  ray-tracing  techniques  must  be  satisfied. 

Part  I  of  this  report  discusses  the  theoretical  aspects  of  the  NISSM  II  com¬ 
puter  program.  Part  II,  a  user's  manual,  contains  several  applications  of  that 
program. 
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To  evaluate  the  ray-tracing  integrals  numerically,  the  velocity  of  sound 
must  be  expressed  as  a  known  function  of  depth.  In  practice,  discrete  data  points 
are  asod  to  generate  an  interpolation  function  from  which  the  velocity  at  mi 
arbitrary  depth  may  be  computed.  This  creates  two  possible  approaches;  (1)  the 
resulting  ray-tracing  equations  con  be  integrated  inclosed  form  or  (2)  the  inter¬ 
polation  function  can  be  used  in  ."exjunction  with  a  numerical- integration  formula. 


The  best-known  variation  of  th*  first  approach  is  the  piecewise  constant- 
gradient  technique,  where  the  sound  velocity  inproxi  mated  by  linear  segments 

as  illustrated  in  figure  2.  The  corresponding  ^ay  diagram,  figure  3,  consists  of 
circular  arcs.  Not.e  the  ray  that  becotueo  ho.  izontal  at  the  velocity  maximum.  It 
is  called  a  splk-iaybeam  since  it  can  be  refracted  upward  as  well  as  downward. 
When  the  boundary-reflected  wys  are  neglected,  the  constant-gradient  technique 
predicts  that  rr  s.rergy  will  eui, .  i.-e  shadow  zone.  This  is  pessimistic,  to  say 
the  least. 


On  tht  oJber  hand,  examples  can  be  constructed  for  which  the  constant- 
gradient  techr’que  is  optimis*  ■  r.  Hern  ring  back  to  equation  (5),  one  sees  that 
the  pressure  -mpdtude,  }  p  J  bfifurxjs>  irfinite  wherever  dx/dCy  vanishes. 
Such  singular  points  a.n  called  caust.  s.  However,  the  constant-gradient  tech¬ 
nique  may  also  predict  false  caustics*'  sand,  thereby,  regions  of  infinite  intensity. 


VELOCITY  (km/») 

1,4  1,5  1.6 


Figure  2.  Piecewise  Constant-Gradient 
Approx;  .nation  of  a  Sound- 
Velocity  Profile 
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RANGE  (km) 


Figure  3.  Ray  Diagram  Plotted  Using  Piecewise 
Constant-Gradient  Techniques 


Another  variation  of  the  first  approach  has  been  implemented  in  NISSM  II; 
that  is,  the  sound  velocity  is  approximated  by  a  continuous  function  of  depth, 
which  has  a  continuous  gradient  and  allows  one  to  integrate  equations  (2),  (3), 
and  (4)  in  terms  of  elementary  transcendental  functions.  Thus  false  caustics 
associated  with  the  constant-gradient  technique  have  been  eliminated,  and  further 
analytical  investigation  is  possible. 

The  velocity  interpolation  function  used  in  NISSM  II  is  given  by 


c(Z)  » 


where 

Az  =  z  -  Zq)  (9) 


+  A  z 


Sq  ' >  +  Az  sf  ' 

(i)  2 

(1  +  Azgp 


V-l/2 


(8) 


in  the  interval 


<  z  <  z 


(i+1) 


(10) 
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Parameters  v^,  g^,  gj*\  and  g^  are  chosen  so  that 

•  c(z)  is  continuous  at  the  velocity  breakpoints, 
cjc 

•  (z)  is  continuous  at  the  velocity  breakpoints,  and 

•  the  maximum  deviation  between  c(z)  and  the  input  data  is  not  to  exceed 
0. 2  m/s. 

The  continuous-gradient  curve-fitting  technique  is  derived  in  detail  in  ref¬ 
erence  3  and  is  summarized  in  appendix  A.  Similarly,  the  range,  time,  and 
range-derivative  integrals  are  evaluated  in  detail  in  reference  4  and  the  results 
are  summarized  in  appendix  B. 


TEMPERATURE-SALINITY-DEPTH  TO  VELOCITY  CONVERSION 

NISSM  II  allows  the  table  of  sound  velocity  to  be  entered  directly,  or  indirectly, 
in  terms  of  temperature,  salinity,  and  depth.  In  the  latter  case,  equation  (7)  in 
reference  5  is  used  to  compute  velocity  from  the  temperature-salinity-depth 
data*;  that  is,  ' 

c  =  c  +  c  +  c.  +  c  +  c ,  ,  (11) 

o  a  b  c  d 


where 


Cq  =  1493.0  +  3{T-  10)  -  6xl0"3(T-10)2  -  4xlO~2(T-  18)2  +  1.2(S-35) 

-  10"2(T  - 18)  (S  -  35)  +  Z/61  ; 

c  =  +  10"1  r2  +  2xl0"4  f2(T-18)2  +  10"1  fv>/90  ;  (13) 

cl 

cb  =  +  2.  6  x  10~4  T  (T  -  5)  (T  -  25)  ;  (14) 

cc  *  -  10"3  f2(f-4)  (f-8)  ;  (15) 

cd  =  1.  5x  10~3  (S-35)2(l-  f )  +  3  x  10“6  T2  (T  -  30)  (S  -  35)  ;  (16) 


*Leroy5  demonstrates  that  his  equations  fit  Wilson's  data0  better  than 
Wilson's  second  formula. 
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c  is  the  velocity  in  meters  per  second, 

Z  is  the  depth  in  meters, 
r  is  the  depth  in  kilometers, 

S  is  die  salinity  in  parts  per  thousand, 

T  is  the  temperature  in  degrees  centigrade,  and 
is  the  latitude  in  degrees. 

The  last  term  in  equation  (11),  c^,  i3  a  corrective  term  for  low  salinities 
and  should  not  be  used  if  S  is  greater  than  30  %0. 


CORRECTION  FOR  EARTH'S  CURVATURE 

The  ray-tracing  integrals  (2),  (3),  and  (4)  were  derived  using  a  rectilinear 
coordinate  system.  Since  it  is  customary  to  express  results  in  terms  of  range 
along  the  ocean  surface  and  depth  below  sea  level,  an  appreciable  error  may 
occur  in  the  ray  trajectories  unless  the  differences  in  coordinate  systems  are 
considered.  For  example,  in  the  rectilinear  coordinate  system  of  figure  4,  a 
ray  in  an  isovelocity  medium  would  be  a  straight  line,  and  the  ocean  surface 
nearly  circular.  However,  if  one  uses  the  range-depth  coordinate  system  of 
figure  5,  the  ocean  surface  becomes  horizontal  and,  to  first  order,  the  ray  is 
an  arc  of  a  circle. 


Figure  4.  Ray  in  an  Isovelocity  Medium 
Plotted  Using  Rectilinear  Coordinates 
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OCEAN  SURFACE 


pv  <JU ,  r.  R  ty  in  an  Isovelocity  Medium 
t -jotted  U'siB?  Range-Depth  Coordinates 


in  NISSM  E,  correction  lor  the  earth's  curvature  is  approximated  by  mo*- 
fying  the  velocity  of  sound  according  to  the  equation 

cmod  “  C(1  +  »  *17) 


where 

c  .  is  the  modified  velocity  of  sound, 
mod 

c  is  the  unmodified  velocity  of  sound, 

z  is  the  depth  of  the  data  point,  and 

R  is  the  mean  radius  of  the  earth  and  equal  to  6370  km. 

Pekerts7  in  bis  discussion  of  e  m"h 

the  errors  introduced  are  small,  providing  that  the  horizontal 

half  the  earth's  radius. 

PROPAGATION  LOSS  ALONG  A  RAY 

The  relative  acoustic  pressure  along  a  ray  is  attained  by  muitiplylug  the 
relative  gemnetrimiTacoustic  pressure  in  equation  ,6,  by  surface,  bottom,  and 
absorption-loss  factors.  Thus 

,  fsur/  fbot  -iwt  +  iwT  +  i*  +  «s  ,  (18) 

P  =  7sp(nsur)  (7botj 
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where 


P  =  -20  log  |  p |  is  the  propagation  loss  in  decibels, 

=  -20  log  j»  is  the  geometrical  spreading  loss  in  decibels, 

N  =  -20  log  n  is  the  surface  loss  in  decibels  per  bounce, 
sur  'sur 

Nfaot  -  “20  log  is  the  bottom  loss  in  decibels  per  bounce, 

A  =  -20  a/ln(10)  is  the  absorption  coefficient  in  decibels  per  kilometer, 
n  is  the  number  of  surface  bounces, 
n^Qt  is  the  number  of  bottom  bounces, 
t  is  the  time  in  seconds. 


T  is  the  travel  time  in  seconds, 


4>  is  the  accumulated  phase  shift  in  radians, 
s  is  the  arc  length  in  kilometers,  and 


w  is  the  frequency  in  radians  per  second. 


SURFACE  LOSS 

Surface  loss  may  be  entered  as  a  table  of  surface  loss  per  bounce,  Nsur, 
versus  the  angle  of  incidence,  0sur.  Here  0sur  is  the  angle  between  the  ray 
and  the  ocean  surface  at  the  point  of  intersection  as  shown  in  figure  6.  The  loss, 
entered  in  decibels  perbounce,  is  converted  to  a  pressure  ratio  according  to 

-N  ur/2°  “Nsur  [ln(10)/20] 

W=1°  =e  (19> 

If  the  surface  loss  is  not  specified,  ijsur  is  set  equal  to  unity  (Ngur  =  0) 
by  default.  The  reader  is  referred  to  the  user's  manual  (Part  II)  for  additional 
details. 

It  must  be  pointed  out  that  the  surface  loss  discussed  above  only  applies  to 
the  propagation-loss  equation  (18).  A  different  surface  loss,  used  in  the  Acous¬ 
tic  Meteorological  and  Oceanographic  Survey  (AMOS)  surface-duct  model,  will 
be  discussed  shortly. 
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OCEAN  SURFACE 


Figure  6.  Intersection  of  a  Ray  with 
the  Ocean  Surface 


BOTTOM  LOSS 


Bottom  loss  may  be  entered  as  a  table  of  bottom  loss  per  bounce,  N^, 
versus  angle  of  incidence,  fl]^,  or  computed  using  either  the  Naval  Undersea 
Center's  (NUC)  bottom- loss  model  in  reference  8  or  the  one  developed  by 
Podeszwa  (reference  9)  at  the  Naval  Underwater  Systems  Center  (NUSC).  All 
three  models  are  functions  of  0^,  asdefinedin  figure  7;  in  each  case,  bottom 
loss  is  converted  io  a  pressure  ratio  according  to 

-Nbot/2°  -Nbot[ln(10)/20] 

"bot =  10  =  6  •  <20> 


■INCIDENT  RAY 


REFLECTED  RAY- 


-OCEAN  BOTTOM 


Figure  7.  Intersection  ox  a  Ray  with 
the  Ocean  Bottom 
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The  NISSM  II  variation  of  the  NUC  bottom-loss  model  is 

Nbot=  [3.7  +  17.5(P-0.27)]  f1/3  .  jtanh  [(®botP/0. 24)1-  5/p] 
+  [(l-P/0.27)/!2.5]  .  (t?bot/90)2J  , 


(21) 


where 

#bot  is  the  angle  of  incidence  in  degrees, 
f  is  the  frequency  in  kilohertz,  and 
P  is  the  porosity. 

This  model  is  evaluated  at  various  frequencies  and  porosities  in  figures  8 
through  11. 

Podeszwa's  bottom-loss  model9  is  an  approximation  of  the  Marine  Geo¬ 
physical  Survey  (MGS)  Acoustic  Provinces  1  through  5.  The  losses  for  that 
model  are  given  by 


L  =2. 2435  ln(0. 1260  6,  .  +  1.496) 
bot  '  bot 

for  Province  1, 

(22) 

lb  =  3.4315  ln(0. 1056  «b  +  2.  842) 

for  Province  2, 

(23) 

tbofc  =  2.4910  ln(0,  8864  $  +  10.  526) 

for  Province  3, 

(24) 

f,  =  2.  8377  ln(l.  8754  6 .  .  +  15.  685) 
bot  bot 

for  Province  4,  and 

(25) 

fbofc  =  2. 4036  ln(20.  576  ebot  +  82. 440) 

for  Province  5 

(26) 

and  are  displayed  in  figure  12.  The  user  of  NISSM  II  is  allowed  to  input  non¬ 
integral  MGS  provinces  between  1  and  5.  A  linear  interpolation  is  performed  in 
determining  coefficients  c^,  C2»  and  c3  for  the  bottom-loss  model, 


Nbot  =  C1  k  (C2  9b*  "  c3>  • 


(27) 


If  the  bottom  loss  is  not  specified,  0bo£  is  set  to  unity  (Nbot  =  0)  by 
default. 
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ABSORPTION  COEFFICIENT 


The  absorption  coefficient  may  be  inputted  as  a  table  of  absorption  per  unit 
distance,  A,  versus  frequency,  f,  or  computed  using  either  Thorp's  equation, 


r  o.i  40 

=  f2  -  + - . 

1  +  f2  4100  +  f2_ 


or  the  Hall-Watson  equation, 11 


A  =  <1. 776  f1* 5 


where 


0. 65053  fT  0. 026847 

+  f 5  - -  + -  /[32.768  +  f3]  ,  (29) 

.  «*  +  4  h  \)j 

[30T  +  102~| 

|_5T  +  2297J 


fT  =  21.  9  x  10 


In  equations  (28),  (29),  and  (30) 

A  is  the  absorption  coefficient  in  decibels  per  kiloyard, 
f  is  the  frequency  in  kilohertz, 
is  the  relaxation  frequency  in  kilohertz,  and 
T  is  the  average  ocean  temperature  in  degrees  Fahrenheit. 


Since  the  units  of  a  in  equation  (18)  are  expressed  in  inverse  kilometers,  then, 
if  either  equation  (28)  or  (29)  is  used  to  compute  A, 


a  =  -  A  [ln(10)/20]  /0.  9144 


When  the  absorption  coefficient  is  not  specified,  A  is  computed  according 
to  Thorp's  equation. 10  (See  figure  13  for  a  comparison  of  the  absorption  coef¬ 
ficient  for  the  two  models. ) 
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FREQUENCY  (kHz) 

Figure  13.  Absorption  Coefficient  versus  Frequency  for  the  Two  Models 

EIGENRAYS 

The  propagation  equation  {(18))  describes  the  variation  in  pressure  along  a 
particular  ray.  Since  a  denumerable  set  of  rays  may  pass  through  a  single  point, 
the  total  pressure  is  the  result  of  the  individual  contributors.  Rays  that  originate 
at  one  given  point,  henceforth  called  the  sonar,  and  terminate  at  another  given 
point,  to  be  referred  to  as  the  target,  are  known  as  eigenrays.  The  problem  is 
to  determine  the  most  important  eigenrays  as  efficiertly  and  as  accurately  as 
possible. 

In  NISSM  II,  eigenrays  are  found  by  tracing  a  preselected  fan  of  test  rays 
to  th<.  target  depth.  When  two  adjacent  rays  of  the  same  type  bracket  a  target 
range,  a  cubic  interpolation  is  performed  to  determine  the  eigenray.  The  two 
steps  for  determining  the  most  important  eigenrays  are  as  follows: 

a.  A  table  of  128  vertex  velocities  (the  cv  of  equation  (1))  should  be  con¬ 
structed.  The  table,  which  begins  with  the  sonar  velocity,  includes  (1)  velocities 
at  breakpoints  in  the  velocity-depth  profile,  at  target  depths,  and  at  the  ocean 
surface  -  ad  bottom  (providing  they  exceed  the  sonar  velocity)  and  (2)  approx¬ 
imately  32  entries  that  also  correspond  to  vertexing  rays.  The  remaining  veloc¬ 
ities  correspond  to  nonvertexing  rays,  the  steepest  of  which  has  an  initial 
inclination  angle  slightly  greater  than  89.  6  deg. 
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b.  One  ray  segment  from  the  sonar  depth  to  the  surface,  another  from  the 
sonar  depth  to  the  bottom,  and  a  third  from  the  sonar  depth  to  the  target  depth 
should  be  traced  for  each  vertex  velocity  in  the  table.  Shallow  ray  segments  are 
terminated  at  their  vertexing  depths.  Each  r;inge,  travel  time,  and  range- 
derivative  increment  is  stored  for  future  reference. 

After  initializing  various  index  counters,  two  test  rays  having  adjacent  ver¬ 
tex  velocities  and  the  same  ray  type  are  generated  by  combining  the  previously 
computed  ray  increments.  The  first  four  ray  types  are  shown  in  figure  14.  Sub¬ 
sequent  types  are  obtained  by  adding  complete  cycles  to  the  basic  four  types.  If 
the  final  ranges  of  the  test  rays  bracket  the  target  range,  a  cubic  interpolation 
is  performed  to  determine  the  eigenray.  If  they  are  not  bracketed,  another  ray 
type  is  considered.  Eventually,  the  end  points  of  adjacent  rays  will  extend  well 
beyond  the  target  range.  Adding  additional  cycles  will  only  move  the  end  points 
farther  to  the  right.  Therefore,  the  ray-type  counter  is  initialized  and  the  ver¬ 
tex  velocity  index  is  incremented. 

After  eight  eigenrays*  have  been  determined,  they  are  sorted  according  to  the 
sonar-inclination  angle  and  entered  into  the  propagation-loss,  the  reverberation, 
or  the  target-echo  routine.  Following  the  processing  of  the  eigenrays,  both  the 
vertex-velocity  and  ray-type  indexes  are  initialized  and  the  target  range  is  in¬ 
cremented. 

Since  the  ray  segments  are  computed  only  once,  the  method  outlined  above 
is  extremely  efficient.  Usually,  because  shallower  rays  tend  to  sustain  fewer 
losses,  the  most  important  eigenrays  are  determined. 


Z 


Figure  14.  Principal  Ray  Types 

*If  a  sound  channel  exists,  up  to  32  eigenrays  may  be  determined. 
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In  NISSM  II,  limiting  rays  (if  they  exist)  are  traced  for  the  first  32  ray 
types  (see  figure  14).  If  the  target  lies  in  the  shadow  of  a  particular  type,  the 
pressure  of  the  corresponding  imaginary  ray  is  computed  according  to  equations 
(32)  and  (33).  *  The  travel  time  of  the  imaginary  ray  is  obtained  by  adding 
csur(x"  xb)  *°  travel  time  of  the  limiting  ray.  Only  three  imaginary  rays 
are  allowed  for  any  particular  target.  Imaginary  rays,  once  determined,  are 
treated  the  same  as  the  other  rays. 

♦For  bottom  shadows,  cgur  and  ggur  are  replaced  by  the  velocity  and 
velocity  gradient  at  the  bottom. 
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Figure  17.  Propagation  into 
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The  procedure  described  above  is  not  entirely  correct.  In  particular,  the 
value  5.  93  in  equation  (33)  was  derived  by  solving  a  free-surface  boundary- value 
problem  that  does  not  apply  to  the  ocean  bottom.  However,  it  is  hypothesized 
that  the  correct  value  should  be  in  the  neighborhood  of  5.  93;  therefore,  equa¬ 
tion  (32)  is  still  meaningful.  The  correct  value  can  easily  be  inserted  when  it 
becomes  available. 


UNIFORM  ASYMPTOTIC  EXPANSION  AT  A  CAUSTIC 

Shadow  zones  are  not  always  caused  by  boundary  surfaces.  Under  suitable 
conditions,  rays  of  a  particular  type  may  be  refracted  away  from  an  ocean  region. 
The  envelope  of  these  rays,  which  is  called  a  caustic,  also  forms  the  boundary 
of  a  shadow  zone. 

In  practice,  caustic  curves  assume  many  shapes;  they  may  even  reduce  to 
a  point  in  a  perfectly  focusing  situation.  However,  we  shall  restrict  this  dis¬ 
cussion  to  that  of  a  smooth  caustic.  The  velocity  profile  in  figure  18  will  produce 
the  smooth  caustic  shown  in  the  ray  diagram  in  figure  19.  Note  that  two  rays 
pass  through  each  point  of  the  illuminated  region  in  the  vicinity  of  the  caustic. 
As  the  caustic  is  approached,  dx/3cv  vanishes  for  both  rays  and  the  geometrical 
spreading-loss  factor,  given  by  equation  (5),  increases  without  bound. 

A  more  accurate  description  of  the  pressure  in  the  vicinity  of  a  caustic  is 
given  by  Ludwig. 13  We  have  incorporated  his  equations  into  the  following  pro¬ 
cedures:  let  T+,  T_,  |p+|,  |p_|  be  the  travel  times  and  pressure  amplitudes 
of  both  the  ray  that  has  touched  the  caustic  and  the  ray  that  has  not  touched  the  caus¬ 
tic,  respectively.  Next  define  and  p  using,  respectively,  equations  (34)  and  (35): 

*  =  (2»f)  (T+  -  TJ/2  (34) 

and 


It  can  be  shown  that  vanishes  on  the  caustic  and  is  positive  in  the  illuminated 
region.  Then,  if 
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a.  2ir  <  i  ,  set 


and 


|p-mod| 

1 

i 

O, 

II 

(36) 

J  P+mod| 

l-N  5 

(37) 

b.  0.  01  <  t  <  1r  ,  set 

1  p-modl 

;  =  Vif[Ai(-p)  (|p+|  +  |p_|)2n/p 

(38) 

+  Ai'(-p)(jp+|  -  | P_|)2/Vp‘J1^2  and 

|  P+mod 

=  0  ; 

(39) 

c.  0  <  ^  <  0.  01,  set 

.  -.1/2 

lp-mod 

|  =  >/F|Ai(-p)  (|p+|  +  |p_|r  Vp  J  and 

(40) 

|pi-mod 

|=0. 

(41) 

Inequations  (38)  and  (40),  Ai(t)  is  the  Airy  function  of  t,  Ai'(f)  is  its 
derivative  and  |P+m(Xj|  and  |P  mo(jj  are  modified  pressure  amplitudes. 

Since  Ludwig's  caustic  expansion^  is  uniform,  in  theory,  one  can  extend 
the  bounds  of  the  inequality  0. 01  <  t  <  2?r  indefinitely,  thereby  eliminating 
cases  (a)  and  (c).  However,  equation  (38)  is  a  lengthier  computation  than  (36) 
and  becomes  numerically  ill-conditioned  as  ’/'—*■  0.  Also,  equation  (38)  can  be 
continued  analytically  into  the  shadow  zone.  In  that  event,  p  is  negati  *e  and 
I  P-mod  |  would  decay  exponentially  with  distance  away  from  the  caustic.  Since 
the  decay  is  extremely  rapid  at  active  frequencies,  the  logic  to  perform  the 
analytic  continuation  does  not  appear  in  the  present  version  of  the  computer  pro¬ 
gram.  This  causes  the  erroneous  discontinuity  in  pressure  indicated  in  figure  20. 


PHASE  SHIFTS 

A  ray  that  either  reflects  off  a  boundary  surface  or  touches  a  caustic  curve 
undergoes  a  corresponding  change  in  phase.  The  accumulation  of  these  abrupt 
phase  shifts  is  signified  by  4>  in  equation  (18).  In  NISSM  II,  <J>  is  initially  set 
to  zero  at  the  sonar  and  decreased  by  r  radians  at  surface  intersections  and 
by  ir/2  radians  at  caustics.  Phase  shifts  at  bottom  intersections  are  set  equal 
to  zero  for  lack  of  data. 
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SURFACE-DUCT  PROPAGATION 

When  the  sound  velocity  initially  increases  and  then  decreases  with  depth 
from  the  ocean  surface,  the  region  of  increasing  velocity  is  called  a  surface 
duct.  Because  ray  theory  does  not  make  allowances  for  leakage  from  the  layer 
into  the  underlying  medium,  and  since  the  ray  geometry  is  more  complicated 
than  that  of  .  a  smooth  caustic,  the  application  of  ray  theory  to  surface-duct 
propagation  is  questionable. 

In  order  to  reduce  the  possibility  of  extremely  erroneous  predictions,  propa¬ 
gation  los3esfor  surface-duct  eigenrays  are  computed  according  to  the  modified 
AMOS  model  in  reference  14  and  ray  theory.  (The  equations  from  that  publi¬ 
cation  are  repeated  in  appendix  C  for  the  reader's  convenience.)  Since  the 
above  paragraph  implies  that  ray-tracing  predictions  are  optimistic,  the  lower 
pressure  (greater  propagation  loss)  is  taken  to  be  the  correct  value.  Only  two 
surface-duct  eigenrays  per  target  are  determined,  and  double  the  AMOS  propa¬ 
gation  loss  factor  (if  used)  is  associated  with  each. 

PROPAGATION  LOSS  VERSUS  RANGE  AT  A  CONSTANT  DEPTH 

One  of  the  most  useful  outputs  provided  by  NISSM  H  is  a  plot  of  effective 
propagation  loss  versus  range  at  constant  depth.  The  value  at  a  particular  range 
is  obtained  by  multiplying  the  relative  acoustic  pressure  of  each  eigenray,  p(*), 
by  factors  corresponding  to  the  reference  level,  dson  or  df.ar,  and  the  respective 
array  responses  j?son-  and  jjtar  The  results  are  then  combined  using  coherent- 
or  random-phase  addition.  Both  sonar  and  target  responses  are  entered  as  tables 
of  response  in  decibels  versus  inclination  angle.  The  reference  level  is  called 
the  sonar  or  target  level,  depending  on  whether  the  system  is  active  or  passive, 
respectively.  Thus  for  coherent- phase  addition,  * 

Xmax 

peff  ~  S  ^son  ^son  ^son^  ^tar  ^tar ^  dfctr  * 
i  =  l 

and  for  random-phase  addition, 
i 

max 

|peff|  “  S  I ^3on  ^son  (®son^  ^tar  ^tar^  ^tar|  »  (43) 

i  =  1 


*In  NISSM  II,  results  for  coherent-phase  addition  are  indicated  by  broken 
lines  and  propagation  losses  computed  using  random-phase  addition  are  indicated 
by  solid  lines. 
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where 


Pg^  =  -20  log  |  Peff  j  is  &e  effective  propagation  loss  in  decibels, 

D  =  +20  log  d  is  the  reference  level  in  decibels  for  an  active  situation, 

son  °  son 

D^ar  =  +20  log  d  is  the  reference  level  in  decibels  for  a  passive  situation, 

p(*)  =  -20  log  j  p^j  is  the  propagation  loss  of  the  i-th  eigenray. 


N  =  -20  log  »  is  the  sonar  response  in  decibels, 
son  &  'son  * 

^ar  =  "20  log  ij  is  the  target  response  in  decibels, 


d^’  is  the  inclination  angle  at  the  sonar  of  the  i-th  eigenray, 
son 


»<i) 

SOI 

5  tar  is  inc^nation  aJ3gle  at  the  target  of  the  i-th  eigenray,  and 

i  is  the  number  of  eigenrays  joining  the  sonar  to  the  target. 

max 

The  default  value  for  d  d  w  and  ,  Is  uuity  (D  D  N 
aud  N(ar=0). 

REVERBERATION  THEORY* 

Consider  an  acoustic  signal  that  originates  at  a  point  source  at  a  reference 
time  equal  to  zero  and  is  transmitted  through  the  ocean.  A  portion  of  the  signal 
is  backscattered  toward  the  source  as  the  signal  encounters  scatterers  on  the 
ocean  bottom  or  the  surface  or  in  the  ocean  volume.  When  the  second  order  scat¬ 
tering  of  sound  is  neglected,  a  closed  ray  path  from  a  source  to  a  scatterer  and 
back  to  the  source  (figure  21)  can  be  constructed  from  an  incident  ray  (a  ray 
from  the  source  at  point  "a"  to  the  scatterer  at  point  "b")  and  a  backscattered 
ray  (a  ray  from  the  scatterer  back  to  the  source).  Let  the  incident  ray  enter  the 
water  at  time  tQ  and  have  travel  time  tp  and  let  the  backscattered  ray  have 
travel  time  t2-  Then  the  closed  ray  path  will  have  round-trip  travel  time, 

(44) 


t  =  ti  +  t2  , 


and  will  return  to  the  source  at  time 


T  =  t  +  t  . 
o 


(45) 


*The  text  under  this  heading  was  excerpted  from  reference  lr 
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Figure  21.  Closed  Ray  Path 

The  acoustic  pressure  at  the  end  of  the  closed  ray  path  is  given  by 

prev  -  dtra  ^tra  ^a.,  1^  P1  P2  ^rec  ^ 9 a,  2^  v  *  ^ 

where 

D^ra  =  +20  log  dtra  is  the  reference  level  in  decibels  of  the  transmitting 

array, 

P1  =  -20  log  p1 1  is  the  propagation  loss  in  decibels  of  the  incident  ray, 

Pg  =  -20  log  p2  is  the  propagation  loss  in  decibels  of  the  backscat- 

tered  ray, 

Ntra  =  -20  log  ij  is  the  transmitter  response  in  decibels, 

N  =  -20  log  n  is  the  receiver  response  in  decibels,  and 

rsc  r6c 

v  is  a  backs cattering  strength  expressing  the  ratio  of  reflected  pres¬ 
sure  to  incident  pressure  per  scatterer. 
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If  the  acoustic  signal  has  pulse  length  r,  it  follows  that 

0  <  tQ  <  r  (47) 

and  that  the  closed  ray  paths  contributing  to  the  reverberation  pressure  at  time 
T  are  those  with  travel  time 


T  -  r  <  t  <  T  . 


The  corresponding  scattered  will  be  contained  in  region  R.  Let  this  region 
be  partitioned  into  numerous  subregions,  AR^,  in  each  of  which  pjp,  pjjp, 
and  are  representative  values,  and  let  there  be  n(')  scat- 
terers  per  unit  region.  Then,  if  the  resultant  reverberation  pressure  is  the 


random-phase  addition  of  the  individual  contributors, 


d  l2  =  yL  /«<*>  )  c(i)o(i)  «  <a{i)  ) 

Prev|  “  (^tra^tra  '®a,  l'  P1  p2  ^rec  '^a, 2} 


fSl)  AR(1)  ,  (49) 


where 


/>  =  (»(,))2 


is  the  backscattering  strength.  In  the  limit  as  AR''  —  0, 


|prev|2“  /Kra’’tra(9a,l,!>lp2’’rec(9a>2)|2,'dR  ' 


Thus  the  reverberation  problem  reduces  to  integrating  equation  (51).  This 
is  accomplished  by  approximating  equation  (51)  with  a  summation  similar 
to  equation  (49).  Since  NISSM  II  usually  computes  eight  paths  (see  footnote  on 
p.  20)  to  each  representative  scatterer,  there  may  be  as  many  as  8^  -  8  =  56 
differently  oriented  closed  paths  per  scatterer. 
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VOLUME  REVERBERATION 


For  volume  reverberation,  equation  (51)  is  replaced  by  the  double  summation, 

,2 


vol 


=  L 


m 


<W,traWa,l)pl  P2  nrec^a,2) 


(52) 


.  M®,  x,l>  Az®A*  , 

vol  * 


where 

A<f>  is  the  horizontal  be  am  width  in  radians, 
x^>  is  the  range  in  kilometers  of  the  (i,  j)-th  scatterer, 
is  the  depth  in  kilometers' of  the  (i,  j)-th  scatterer, 
and 

Uvol  =  10  log  Mvol  is  volume-scattering  strength  per  cubic  kilometer 
in  decibels. 

The  horizontal  beamwidth,  A4>,  equidistant  ranges,  x(i>,  and  stable  of 
volume-scattering  strengths,  versus  depth,  z<i),  are  input  data.  They 

are  converted  by  NIS3M II  to  the  appropriate  units.  Then  zO),  ••  0),,  and  Az(J) 
are  computed  according  to 


and 


<£«♦»>  ♦J®)  . 

(53) 

(54) 

A4®»Z«+1>-Z®  , 

(55) 

respectively.  At  the  ranges  of  interest,  the  quantity  x^Axtf*  J)az0)a$  is  a 
good  approximation  of  the  volume  of  the  insonifled  region  bounded  by  zG),  z0+1) 
and  surfaces  of  constant  round-trip  travel  time  (see  figure  22). 
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Figure  22.  Insonii.ei  Region 


Initially,  the  range  increment  Ax^’^  is  set  to 
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(56) 


and  the  round-trip  travel  times,  t^’^~,  at  (x^  +  Ax/2^’^,  z^)  are  aonroxi- 
mated  by 


tfcjte.tfc#±  / ‘  »  \  tafcjv,. 


“  CM  c#’»J 

v,l  °v,2, 


(57) 


Kv 


Here  and  are  the  vertex  velocities  of  the  incident  and  backs  cat- 

v,  1  v,  2 


tered  rays,  respectively.  Then,  if 

T-r<tM“  <  t(i,3>  <  t(i,i)%  T  t 


(58) 


(i  ^)4* 

the  region  whose  sides  are  surfaces  of  constant  t  ’  is  contained  within  the 
insonified  region  and  Ax^*  3)  is  left  unchanged, 
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However,  if  the  inequality  (58)  does  not  hold,  modified  values  of  t(M)±  and 
Ax(i»  J)  are  determined  according  to 


•S£--t«  ‘(i'ih.T-r]  , 


tmi!+*[mln  t<1,i)+'  T]  • 


4,0. «,  /J_t_L 

mod  l  (i,j)  (i,j) 

\v,l  Cv,  2 


nr)  •  <61> 


Each  term  in  the  volume-reverberation  summation  ((52))  has  now  been  found. 


If  the  volume-scattering  strength  table  is  omitted,  then  the  volume-rever¬ 
beration  computation  is  also  omitted. 


SURFACE  REVERBERATION 


Surface  reverberation  is  computed  according  to  summation  (49),  with 


AR(i)  =  M(i)  *(i)  Ax{l)A* 
sur 


(i) 


10  log  =  3. 3/3  logl-  —  — — — 2-  )  -42,4  log/3  +  10  log  m  ,  (63) 

SUr  \  30  30  /  SUr 


where 


®sur  1  *s  sur^ace"6razinS  angle,  in  degrees,  of  the  incident  ray  to  the 
’  i-th  surface  scatterer, 


6 w  2  *s  sur^ace"Sraz^nS  angle,  in  degrees,  of  the  backscattered  ray 
*  from  the  i-th  surface  scatterer,  and 


U  =10  log  m  is  a  surface-scattering  constant  per  square  kilometer, 
sar  SUr  in  decibels,  and  /3=  358  (vf1/3)-0*  53, 


where 


_ -- . . .  — - i  ********* 'll  S\V- *■>',>>-"<  'tXi'i'i \,',V "i T' ^  ~''* 
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v  is  the  wind  speed  in  knots  and 
f  is  the  frequency  in  hertz. 

Equation  (63)  reduces  to  the  Chapman-Harris  surf ace-r> catering  strength16 

when  0^  ,  equals  0^  _  and  U  is  set  equal  to  2. 6dB.  If  U  is 

sur,  1  sur,  2  —  sur  —  sur 

omitted  from  the  input  data,  then  the  surface  reverberation  is  also  omitted. 
BOTTOM  REVERBERATION 

Bottom  and  surface  reverberations  are  computed  similarly,  except  that  for 
bottom  reverberation 

10  “«  "tot '  10  log  (sin  C,  1  sin  6  hot,  2)  +  10  log  "bot  '  (64) 

where 

—bot  =  10  log  -bot  iS  a  bottom-scattering  constant  per  square  kilometer 
in  decibels. 

When  U^qj.  equals  -27  dB,  equation  (64)  reduces  to  MacKenzie's  bottom¬ 
scattering  strength. 17  If  is  omitted,  then  the  bottom  reverberation  is 

not  computed. 

TOTAL  REVERBERATION 

A  ray  leaving  the  transmitter  with  an  initial  inclination  angle  of  +90  deg  v/ill 
be  reflected  from  the  ocean's  surface  and  bottom  and  will  return  energy  to  the 
receiver.  This  phenomenon  is  observed  in  actual  measurements.  The  pressure 
associated  with  such  eigenrays  is  computed  using  the  propagation-loss  equa¬ 
tion  ((43))  with  subscripts  "son"  (sonar)  and  "tar"  (target)  replaced  by  "tra" 
(transmitter)  and  "rec"  (receiver),  respectively.  The  total  reverberation  is 
assumed  to  be  the  random- phase  addition  of  surface  and  bottom  echoes  and 
volume,  surface,  and  bottom  reverberations. 
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TARGET  ECHO 


The  target  echo  at  a  particular  time  is  determined  using  the  summation 

o  i2=y  d  „  io(i)  i o(i) D(i) »  (o{i>  )2 u{i)  « 

peeho|  ~  tra  ^tra '®a,  1^  P1  P2  \ec '®a,  2'  ^tar  ’ 


where  U.  =  10  log  £  is  the  target  strength  in  decibels. 

—tar  tar 

As  for  reverberation,  the  summation  above  is  taken  over  all  closed  paths, 
with  round-trip  travel  times  between  T  -  t  and  T.  If  U  is  omitted,  then 
the  target  echo  is  not  computed. 


NOISE 


It  is  convenient  to  separate  noise  into  the  two  components: 


•  P  ^  =  20  log  p  ^  is  the  self-noise  spectrum  density  in  decibels,  and 

•  pamj;)=  20  log  Pamb  is  the  ambient-noise  spectrum  density  in  decibels. 


A  single  value  for  self-noise,  Pseif,  may  be  included  in  the  input  data.  If 
not  specified,  pgeif  is  assumed  to  be  zero. 

Ambient  noise  maybe  entered  as  a  table  of  Pamb  versus  frequency  or  taken 
from  an  ambient-noise  curve,  which  is  obtained  by  adding  the  contributions  from 
the  moderate  shipping  and  wind  speed  shown  in  figure  23.  This  graph,  obtained 
by  approximating  figure  7,  5  in  reference  18,  contains  the  Knudsen  spectra  in  the 
500-  to  50, 000-Hz  band.  The  default  value  for  ambient  noise  is  determined  by 
evaluating  the  stored  curve  at  the  particular  frequency  of  interest. 

Since  self-noise  and  ambient  noise  are  assumed  to  be  isotropic,  the  total 
ncise-spectrum  density  at  the  beamformer  output  is 

P»oise  *  10  l0g  [”DI  <PLf  +  pamb>]  <66> 

where  DI  =  -10  log  >?Dj  is  the  directivity  index  in  decibels.  The  default  value 
for  t?di  is  unity  (DI=0). 
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Figure  23.  Ambient  Noise  versus  Frequency 
(reference  18) 


100,000 


TARGET  ECHO  TO  MASKING  BACKGROUND 

For  a  narrowband  process,  the  signal-to-noise  ratio  (s/n),  that  is,  the  target 
echo  to  masking  background,  at  a  particular  time  T  may  be  approximated  by 


s/n  = 


Af  p 
echo  echo 


2  2 
Af  .  p  .  +  Af  p 

noise  noise  rev  rev 


where 


Pecho  =  20  Pecbo  is  tarSet  ec^°  spectrum  density  in  decibels  at 

time  T, 

P  =20  log  p  is  the  reverberation  spectrum  density  in  decibels  at 
rev  rev  _ 

time  T, 

P  .  =  20  log  p  .  is  the  noise- spectrum  density  in  decibels  at  the 

noise  noise  ,  ,  ,  . 

beamformer  output, 


SshjS^&'a,.:  'iv&kaik?','  iliXj 
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PROBABILITY  OF  DETECTION 

Probability  of  detection  is  a  useful  concept  in  evaluating  the  performance  of 
sonar  systems.  The  particular  narrowband  detector  modeled  in  NISSM  II  and  de¬ 
picted  in  figure  24  was  developed  by  Nuttall  and  Magaraci.  i9 


■SQUARE-LOW  ENVELOPE  DETECTOR 
Figure  24.  Narrowband  Detector 


The  input  x(t)  is  assumed  to  be  either  a  stationary  zero-mean  Gaussian 
signal,  s(t),  and  noise,  n{t),  or  noise  alone.  It  is  assumed  that  the  signal  spec¬ 
trum  has  the  same  B  hertz  bandwidth  as  the  narrowband  filter  and  is  centered 
on  it,  and  that  the  noise  is  flat  over  the  frequency  interval.  Samples  of  the 
squared-envelope,  which  are  taken  every  B-1  seconds,  are  accumulated  for  an 
observation  time  of  T  seconds.  The  threshold,  fl  ,  is  fixed. 

Then,  if  PF  is  a  given  probability  of  a  false  alarm,  a  quantity  A  is  deter¬ 
mined  by  inverting  the  expression 


(72) 


where  S/N  is  the  maximum  value  with  respect  to  time  of  s/n  (equation  (67))  at  a 
particular  target.  If  the  probability  of  a  false  alarm  is  omitted  from  the  input 
data,  then  the  probability  of  detection  is  not  computed. 
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REMARKS 

NISSM  n  is  a  combination  of  numerous  theoretical  and  empirical  models 
joined  into  one  highly  sophisticated  computer  program.  It  is  unfortunate  that  the 
developers  of  these  models  chose  to  express  their  results  in  so  many  different 
sets  of  units.  Although  the  computer  execution  time  required  to  convert  data  to 
a  consistent  set  of  units  is  inconsequential,  a  great  deal  of  otherwise  unnecessary 
logic  must  be  added  to  the  program. 

We  chose  to  express  our  equations  in  a  form  that  would  approximate  the 
original  model  as  closely  as  possible  so  that  the  programmer  and  user  could 
easily  compare  the  original  model  with  the  one  presented  here. 

Many  of  the  environmental  routines,  such  as  bottom  loss  and  attenuation, 
can  easily  be  replaced  as  more  accurate  formulas  are  developed  Similarly,  one 
may  wish  to  use  a  more  applicable  detection  model.  We  have  tried  to  anticipate 
these  as  well  as  other  possible  modifications  and  have  designed  NISSM  II 
accordingly. 
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Part  II 

WZfVH  MANUAL 

This  manual  dhscr  ibes  the  input  requirements  of  NISSM  n.  It  is  assumed 
■hat  ‘Me  reader  has  read  Part  t  and  is  familiar  with  the  terminology  defined 

'herein. 


ffARDWARE  AND  SOFTWARE  REQUIREMENTS 

At  present,  NISSM  11  is  operational  on  both  the  NUSC  and  NUC  Univac  1108, 
but  there  are  miner  differences  in  the  NUSC  and  NUC  models: 

•  The  FORTRAN  V  compilers  are  not  compatible  with  respect  to 
ENCODE-DECODE  statements.  The  arguments  for  these  state¬ 
ments  are  in  reverse  order. 

•  The  NUSC  CalComp  plot  routines  AXIS  and  SCALE  nave  one  more 
argument  than  those  of  NUC. 

For  both  models  the  input  is  read  from  SO- column  punched  cards,  and  no 
peripheral  equipment,  save  the  Cal  comp  batch  plot  tape  drive,  is  needed  to  store 
intermediate  results. 

Table  1  lists  the  various  subroutines  that  constitute  the  computer  program, 
but  does  not  include  the  standard  library  software  —  it  was  derided  to  simplify' 
the  logic  or  use  similar  logic  whenever  feasible,  even  at  the  expense  of  com¬ 
putational  efficiency.  The  programmer,  once  he  has  mastered  a  particular 
environmental  subroutine,  should  be  able  to  modify  the  remaining  subroutines 
without  much  difficulty.  Approximately  100,  000g  words  of  core  are  required. 
This  number,  which  includes  standard  library  software,  could  be  reduced  by 
decreasing  the  dimensions  of  several  large  arrays. 
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Table  1.  Subroutines 


Name 

Entry  Point 

Description 

AIRY 

Computes  both  kinds  of  Airy  functions  and  their 
derivatives  for  small,  real  arguments. 

AMB5NT 

Computes  ambient  noise  by  interpolating  in  a 
table  of  noise  versus  frequency  or  by  using 

Knudsen  curves. 

INAMBN 

Reads  a  table  of  ambient  noise  versus  frequency 
and  sets  options. 

AMOS 

Computes  propagation  loss  in  surface  ducts 
according  to  the  modified  AMOS  equations. 

AMRA 

Determines  the  maximum- response  axis  of  the 
various  arrays. 

ATTEN 

Computes  attenuation  by  interpolating  in  a  table  of 
attenuation  versus  frequency  or  by  using  either 
the  Thorp*0  or  the  Waison-Hall*1  equation. 

INATTN 

Reads  a  table  of  attenuation  versus  frequency  and 
sets  options. 

BOTLOS 

Computes  bottom  loss  by  interpolating  in  a  table 
of  bottom  loss  versus  bottom  angle  or  by  using 
either  MGS  bottom-loss  curves  or  the  NUC 
bottom-loss  model. 

INBTLS 

Reads  a  table  of  bottom  loss  versus  bottom 
angle,  the  MGS  province,  or  the  bottom  porosity 
and  sets  options. 

DERF 

Evaluates  the  error  function  in  double  precision. 

DINERF 

Evaluates  the  inverse  error  function  in  double 
precision. 

FACTR 

Converts  input  units  to  NISSM  n  units  according 
to  table  2 . 
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Table  1  (Coat'd).  Subroutines 


Name 

Entry  Point 

Description 

FITVEL 

Reads  the  velocity  or  temperature-salinity  pro¬ 
file  and  returns  the  corresponding  continuous- 
gradient  velocity  parameters. 

GRADE 

Evaluates  the  gradient  of  the  continuous-gradient 
velocity  profile  at  a  particular  depth. 

INPUT 

Reads  input  data  directly  or  transfers  to  an  entry 
point  in  a  subroutine  that  reads  input  data. 

INTERP 

Interpolates  linearly  in  tables. 

LAYER 

Locates  the  continuous-gradient  velocity-profile 
layer  in  which  a  depth  is  situated. 

NISSM2 

Oversees  subroutines  of  the  NISSM  n  program. 

PRBDET 

Determines  the  probability  of  detection  as  a 
function  oi  oignal-to-noise  ratio  for  a  given 
probability  of  false  alarm  and  sample  rate. 

PRNTIG 

Prints  acoustic  eigenrav  data. 

INPRIG 

Initializes  subroutine  PRNTIG. 

PRNTRY 

Prints  coordinates  along  rays. 

INPRRY 

Initializes  PRNTRY. 

PROBF 

Evaluates  the  Gaussian  probability  function. 

RAY 

Integrates  the  ray-tracing  equations  over  one 
continuous-gradient  velocity-profile  layer. 

RAYPLT 

Constructs  ray  diagrams. 

RECVR 

Computes  the  response  of  the  receiver  array  by 
interpolating  in  a  table  of  response  versus 
receiver  angle. 
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Table  1  (Cont'd).  Subroutines 


Name 

Entry  Point 

Description 

INRCVR 

Reads  a  table  of  receiver  response  versus 
receiver  angle  and  sets  options. 

REVERB 

Computes  the  contribution  of  a  scatterer  to 
bottom,  surface,  and  volume  reverberation. 

RTNSCL 

Returns  CalComp  scale  factors. 

INAXIS 

Reads  CalComp  axis  information. 

SHAPE 

Computes  the  pulse  shape  of  the  signal  at  a 
target. 

SONAR 

Computes  the  response  of  the  sonar  array  by 
interpolating  in  a  table  of  response  vers  us  sonar 
angle. 

INSONR 

Reads  a  table  of  sonar  response  versus  sonar 
angle  and  sets  options. 

SURLOS 

Computes  the  surface  loss  by  interpolating  in  a 
table  of  loss  versus  surface  angle. 

INSRLS 

Reads  a  table  of  surface  loss  versus  surface 
angle  and  sets  options. 

TABLE 

Reads  tables. 

TARGET 

Computes  the  response  of  the  target  array  by 
interpolating  in  a  table  of  response  versus 
target  angle. 

_ 

INTARG 

Reads  a  table  of  target  response  versus  target 
angle  and  sets  options. 
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Table  1  (Cont'd).  Subroutines 


Name 


Entry  Point 


Description 


TRANS 


Computes  the  response  o£  the  transmitter  array 
by  interpolating  in  a  table  of  response  versus 
transmitter  angle. 


VELOC 


INTRAN  Reads  a  table  of  transmitter  response  versus 
transmitter  angle  and  set  options. 

Evaluates  the  velocity  of  the  continuous-gradient 
velocity  profile  at  a  particular  depth. 


VELTMP 


Converts  temperature-salinity-depth  10  velocity. 


VOLSCA 


INVLSC 


Computes  the  volume-scattering  strength  by 
interpolating  in  a  table  of  scattering  strength 
versus  depth. 

Reads  a  table  of  volume-scattering  strengths 
versus  depth  and  sets  options. 
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Tsvbl-3  2.  Input  Units  and  Conversion  Factors 
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Dimension 

Input  Unit 

Conversion  Factor 

Code 

Results 

Length 

certimeters 

1.  OOOOOOOOxKT5 

CM 

km 

inches 

2.  S4000000xl0-5 

E* 

feet 

0. 30480000X10-3 

FT 

yards 

0.  9 1440000 xlO-3 

YD 

meters 

1.  OOOOOOOOxKT3 

M 

fathoms 

1.  82880000xl0“3 

F 

kilofeet 

0.30480000 

KFT 

kiioyards 

0.91440000 

KYD 

kilometers 

1.00000000 

KM 

miles 

1.60934400 

MI 

nautical  miles 

1.85200000 

N  MI 

degrees  latitude 

1. 11 120000 xl 0+2 

DEG  L 

Velocity 

centimeters/seco:  : 

1.  OOOOOOOOxKT5 

CM/S 

km/s 

inches/second 

2. 54000000xl0“5 

IN/S 

feet/second 

0. 30480000xKT3 

FT/S 

yards/second 

0.  91440000xKT3 

YD/S 

meters/second 

1. OOOOOOOOxlO-3 

M/S 

fathoms/seccnd 

1.  82880000 xlO-3 

F/S 

kilofeet/second 

0.30480000 

KFT/S 

kiloyards/second 

0.91440000 

KYD/S 

kilometers/second 

l. 00000000 

KM/S 

miles/second 

1.60934400 

MI/S 

knots 

0. 514444 44xl0-3 

KNOTS 

Angle 

degrees 

1. 74532925X10-2 

DEG 

rad 

radians 

1.00000000 

RAD 

t 

Frequency 

radians/second 

0. 15915494xl0-3 

RAD/S 

kHz 

hertz 

1. OOOOOOOOxlO-3 

HZ 

cycles/second 

1. OOOOOOOOxlO-3 

CPS 

kilohertz 

1.00000000 

KHZ 

kilocycles/ second 

1.00000000 

KCPS 

' 

Time 

milliseconds 

1. OOOOOOOOxlO-3 

MS 

f 

seconds 

1.00000900 

SEC 

f 

Temperature 

degrees  centigrade 

1.00000000 

DEG  C 

o 

P 

degrees  Fahrenheit 

5/9(T-32) 

DEG  F 

t 

Less 

decibels 

1.00000000 

DB 

dB 

MGS 

1.00000000 

MGS 

i 

mssami 


-  -Tar« 


t  i 
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Table  2  (Cont'd) .  Input  Units  and  Conversion  Factors 


Dimension 

Input  Unit 

Conversion  Factor 

Code 

Results 

! 

Attenuation 

decibels/centimeter 

1.00000000x10+5 

DB/CM 

dBAm 

decibels/inch 

0.39370079x10+5 

DB/IN 

decibels/feet 

3.28083990x10+3 

DB/FT 

decibels/yard 

1. 09361330x10+3 

DB/YD 

decibels/meter 

1. 00000000x10+3 

DB/M 

decibels/fathom 

0.54680665x10+3 

DB/l? 

aecibeis/kilofeet 

3.28083990 

DB/iCFT 

decibels  Ailoyard 

1.09361330 

DB/KYD 

decibels  Ailometer 

1.00000000 

DB/KM 

deoibels/mile 

0.62137119 

DB/MI 

decibels/nautical  mile 

0.53995680 

DB/N  M 

Dimensionless 

thousand  to  1 

1.  00000000x10+3 

1000/ 

hundred  to  1 

1.  00000000xl0+2 

100/ 

ten  to  1 

1.  OOOOOOOOX10+1 

10/ 

1.00000000 

1/ 

1.00000000 

/I 

parts  per  ten 

1.  00000000 xio-l 

/io 

parts  per  hundred 

1.  00000000x10-2 

/100 

parts  per  thousand 

1.  00000000 X10-3 

/1000 

PARAMETER  CARDS 

NISSM  II  input  data  cards  are  grouped  according  to  format  in  table  3A, 
where  we  see  that  parameter  cards  specify  a  single  number  such  as  the  sonar 
depth.  The  fust  30  columns  of  the  card  contain  the  data-identification  cede, 
the  next  10  list  the  value,  and  the  following  10  indicate  the  units.  For  example, 
the  following  parameter  card  indicates  that  t?-e  sonar  depth  is  20  ft: 


/ 


DNAR  DEPTH 


20.0  FT 


\ 


00000000  OuOOOOOOOftOOOOOOOOOOOOOOOO  0  00  oooooooooooooooo 

i  :  3  <  ,  5  i  8  » iq  n  i»  u  h  » it  ini » :i  ini  ins  w  »  a  a  » in;  am  » aum  as « 4i»  m  h  it «  it  ii  h  m  it  « «  u  «  «  « <» 

in  1 1  m  1 1 1  m  1 1  in  1 1  in  i  \\vy  ^4 1 1 1 1 

J  2  2  2  2J 


In  the  majority  of  cases,  the  unit  may  be  any  alphanumeric  abbreviation 
taken  from  table  2,  and  the  input  parameter  is  simply  multiplied  by  the  corre¬ 
sponding  conversion  factor.  Scattering  strengths  and  intensity  levels,  on  the 
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other  hand,  are  converted  to  pressure  ratios  and  must  be  entered  in  decibels. 

A  large  negative  number  less  than  -300  dB  will  be  converted  to  zero. 

If  a  parameter  card  is  omitted  from  the  input  deck,  the  default  value  is 
assumed.  If  two  or  more  parameter  cards  having  the  same  identification  code 
appear  in  the  input  deck,  the  most  recently  read  value  is  used. 

It  is  necessary  to  discuss  several  of  the  input  parameters  in  somewhat 
more  detail:  The  minimum  and  maximum  ranges  define  the  target-range  coor¬ 
dinates  by  generating  128  equally  spaced  values.  Similarly,  the  minimum  and 
maximum  times  generate  256  equally  spaced  values  for  the  reverberation  and 
signal-to-noise  computations. 

We  see  from  table  3A  that  one  may  input  the  sea  state,  ss,  wave  height, 
zw,  and  wind  speed,  v,  as  three  independent  quantities.  However,  if,  for 
example,  the  sea  state  is  specified,  but  the  wind  speed  and  wave  height  are  not, 
then  v  and  zw  are  computed  according  to 

v  =  5  x  ss  knots 
and 

z  =  0.0026  x  v5/2  ft, 
w 

respectively.  The  last  equation  may  be  found  in  reference  20.  If  not  one  of  the 
three  quantities  is  inputted,  all  are  assumed  to  be  zero  by  default. 

Parameter  cards  must  be  punched  as  follows: 

Columns  1-30  Data-identification  code  left  justified. 

Columns  31-40  Value  in  an  F10. 0  format. 

Columns  41-50  Unit  left  justified. 
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Table  3A.  Input  Parameter  Cards 


Data- Identification  Code 
Averaging  Time 


Bandwidth 


Bottom  Depth 

Bottom-Scattering  Constant 
Closing  Speed 
Directivity  Index 
Frequency 

Horizontal  Beamwidth 
Maximum  Range 
Maximum  Time 
Minimum  Range 
Minimum  Time 
Probability  of  False  Alarm 
Pulse  Length 
Receiver  Level 
Receiver  Tilt  Angle 
Sea  State 


Self-Noise 
Sonar  Depth 
Sonar  Level 
Sonar-Tilt  Angie 
Surface-Scattering  Constant 
Target  Depth 


Default  Value 


Pulse  length 


1  Hz  if  the  pulse  length  has  not  been  read; 
(pulse  length)-*-  if  the  pulse  length  has  been 
read. 

Maximum  depth  in  the  velocity-depth  profile. 
0  (-«  dB) 

0  km/s 
1  (0  dB) 

0  Hz 

2jt  rad  (360  deg) 

0  km 


0.1  km 
0.1  s 


1  <0  dB) 

0  rad 

Function  of  wave  height  or  wind  speed, 
depending  on  which  is  read  in  last. 

0  (-00  dB) 

0  km 

1  (0  dB) 

0  rad 

0  (-00  dB) 

0  km 


«« fry , 
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Table  3A  (Cont’d).  Input  Parameter  Cards 


Data- Identification  Code 

Default  Value 

Target  Level 

1  (0  dB) 

Target  Strength 

0  (-00  dB) 

Target-Tilt  Angle 

0  rad 

Transmitter  Level 

1  (0  dB) 

Transmitter  Tilt  Angle 

0  rad 

Wave  Height 

Function  of  sea  state  or  wind  speed,  depend¬ 
ing  on  which  is  read  in  last. 

Wind  Speed 

Function  of  sea  state  or  wave  height,  depend¬ 
ing  on  which  is  read  in  last. 

The  data  in  table  3B  maybe  entered  as  a  table  of  dependent  variables  versus 
independent  variables  according  to  the  format  below: 


Card  1 

Columns 

1-30 

Data-identification  code  left  justified. 

Card  2 

Columns 

Columns 

1-10 

11-20 

Units  of  the  independent  variable  left  justified. 
Units  of  the  dependent  variable  left  justified. 

Card  i+2 

Columns 

1-10 

i-th  entry  of  the  independent  variable  in  an 
F10.0  format.  (The  entries  must  form  a 
strictly  increasing  sequence. ) 

Columns 

11-20 

i-th  entry  of  the  dependent  variable  in  an 

F10. 0  format. 

Last  card 

Columns 

1-5 

7/8  EOF  (end  of  file) 

When  indicated,  a  compact  format  is  allowed.  The  independent  variable 
entries  are  assumed  to  be  spaced  one  degree  apart,  and  only  the  dependent 
variable  entries  are  read: 

Card  1  Columns  1-30  Data-identification  code  left  justified. 

Columns  31-40  Word  COMPACT  left  justified. 

Card  i  +  l  Columns  1-50  Dependent  variable  entries  in  a  10 F5. 0 

format. 
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Axis  cards  listed  in  table  3C  must  be  punched  as  follows: 

Columns  1-30  Data-identification  code  left  justified. 

Columns  31-40  Minimum  value  to  be  plotted  in  an  F10.  0  format. 
Columns  41-50  Maximum  value  to  be  plotted  in  an  FI 0. 0  format. 
Columns  51-60  Axis  length,  in  inches,  in  an  F10. 0  format. 

Columns  61-70  Axis  unit  left  justified. 

The  default  values  listed  below  pertain  to  the  axis  length  only.  Under  a 
default  condition,  the  minimum  and  maximum  values  are  set  by  the  CalComn 
subroutine  SCALE. 


Table  3C.  input  Axis  Cards 


The  comment  card  allows  the  user  to  print  80-column  comments  on  the 
printed  output  and  CalComp  plots: 

Card  1  Columns  1-30  Word  COMMENT  left  justified. 

Card  i+1  Columns  1-80  80-column  comments. 

Last  card  Columns  1-5  7/8  EOF. 

54 


TR  4527 


The  process  card  signifies  that  ail  environmental  data  have  been  read  and 
that  the  main  program  should  begin  processing  the  information: 


Columns  1-30 
Columns  31-40 

Columns  41-50 


Word  PROCESS  left  justified. 

Word  PRINT,  left  justified,  will  cause  eigenray  data 
to  be  printed. 

Word  TOTAL,  left  justified,  will  prevent  the  reverber¬ 
ation  computation  from  being  repeated;  word  VOLUME, 
left  justified,  will  prevent  the  volume-reverberation 
computation  from  being  repeated. 


The  initial  inclination  angle  card  specifies  which  rays  are  to  be  traced  by 
RAYPLT: 


Columns  1-10 
Columns  11-20 
Columns  21-30 
Columns  31-36 
Columns  37-40 


Columns  41-50 

Columns  51-60 

Columns  61-70 
Columns  71-80 


Initial  angle  of  the  first  ray  in  an  F10. 0  format. 

Initial  angle  of  the  final  ray  in  an  F10. 0  format. 

Initial  angle  increment  in  an  FI  0. 0  format. 

Angle  units,  loft  justified. 

A  print  option  in  an  14  format.  If  blank  or  0,  the 
printing  of  ray  coordinates  will  be  suppressed.  If  i, 
the  final  coordinates  of  every  i-th  ray  segment  will 
be  printed. 

Maximum  number  of  surface  intersections  in  an  110 
format. 

Maximum  number  of  bottom  intersections  in  an  110 
format. 

Maximum  number  of  upper  vertexes  in  an  II 0  format. 
Maximum  number  of  lower  vertexes  in  an  II 0  format. 


The  end  card  signifies  that  the  execution  should  be  terminated: 
Columns  1-30  Word  END  left  justified. 
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INPUT  TABLES 

Environmental  and  system  characteristics  that  ordinarily  can  not  be  de¬ 
scribed  by  a  single  parameter  may  be  entered  as  a  table  of  independent  versus 
dependent  variable.  The  standard  table  format  requires  punching  the  data- 
identification  code  in  columns  1-30  of  the  first  card  in  the  card  set.  Card  2 
contains  the  independent  and  dependent  units  in  columns  1-10  and  11-20,  re¬ 
spectively.  Then  the  table  is  inputted  one  data  point  per  card  according  to  a 
2F10. 0  format.  An  EOF  follows  the  last  data  card  of  the  set.  For  example, 
the  following  is  a  table  of  velocity-profile  data. 


42000.0 4946. 0 


9000.0  4793. 0 


9000.0  4788. 0 


bOOTO 4818.0 


f00. 0  4843.0 


f0070 4871.0 
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Velocity  profiles  may  also  be  expressed  in  terms  of  depth,  temperature,  and 
salinity.  Allowable  temperature  units  can  be  in  degrees  Fahrenheit  or  centi¬ 
grade,  but  salinity  is  usually  in  parts  per  thousand.  The  latitude,  in  degrees, 
appears  in  columns  31-40  of  the  units  card.  Depth,  temperature,  and  saliDity 
data  punched  in  columns  1-10,  11-20,  and  21-30,  respectively,  follow.  For 
example, 
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ioio.o  54.0 


4®670 37.4 


40076  36.9 


4676  70.7 


70.9  3A.23 


7178 38.23 
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The  standard  table  format  permits  inserting,  deleting,  and  correcting 
entries  with  minimum  effort.  Independent  variables  must  forj  i  an  increasing 
sequence,  but  are  otherwise  arbitrary.  A  possible  disadvantage  of  the  standard 
table  format  is  that  the  input  deck  becomes  difficult  to  manipulate  when  several 
large  tables  are  present. 

Thus  we  have  given  the  option  to  enter  certain  environmental  and  system 
characteristics  in  a  fixed  compact  format.  The  word  COMPACT  must  be 
punched  in  columns  31-40  of  thedata-identification  card,  which  is  to  be  followed 
by  the  dependent  variables,  in  decibels,  according  to  a  10F5. 0  format.  If  the 
compact  option  is  used,  181  values  of  response  must  be  entered  for  every  de¬ 
gree  from  -90  to  +90.  Surface  and  bottom  losses,  however,  must  be  entered 
for  every  degree  from  0  to  90;  for  example, 


57 


TR  4527 


Several  environmental  characteristics  may  be  generated  internally  by  the 
computer  program.  For  example,  the  following  card  will  activate  the  MGS  3 
bottom-loss  function: 


BOTTOM  LOSS 


3.0  MGS 


00'  00000'  '00300000D000300000000000000  0  D  0  0000000000? 

I  i  J  (  5  $  I  I  3  1C  11121)1115 16  Wit  IS  :>!!;<!$  26  !1?S  IS  .’0)1  32  33  25  !S  ’6  HU  53  <0111!  J1<US<E<MI<9  93 SI  SiiisX 


I  111  1  1  III  I  III 


111  I  III  111  I  I  1  I 


See  table  3B  for  a  list  of  additional  options. 

CALCOMP  AXIS  CARDS 

The  user  may  either  (1)  set  CalComp  axis  values,  lengths,  and  units  or  (2) 
use  default  values.  If  the  first  option  is  employed,  the  appropriate  axis  card 
must  contain  the  data-identification  code  in  columns  1-30,  the  minimum  and 
maximum  values,  the  axis  length,  and  the  output  units  in  the  four  subsequent 
field  of  10  columns  each.  The  following  axis  card  will  allow  reverberation 
values  from  -90  to  0  dB  to  be  plotted  on  a  graph  whose  vertical  axis  is  9  in. 

/fverberatidm  axis  “  =90  oTo  976  db~7 


00  000000  00030  0  OCOflOOOOOflOOOO  0  00000  0  000080000  000000000? 

I  i  j  3 «  s  >  i  1 9  wi;  i.’ mi  ims '.i  muon  nnunnunnm\  mm»»3?  n  » <a o  u  <s  <unn9  to  si » S3  ss  s(  sm  si  m  h  s; » 

N.1HHH  n m1 1 i ii m 1 1 1 1 mj>" — 

Default  values  for  axis  lengths  are  fisted  in  table  3C;  minimum  and  maxi¬ 
mum  default  values  are  set  by  the  CalComp  subroutine  SCALE,  and  the  default 
units  are  determined  by  other  input  data  cards. 


PROCESS  CARD 

The  process  card  serves  to  separate  oases  and  set  certain  options.  If  the 
word  PRINT  appears  in  columns  31-40,  the  acoustic  eigenrays  to  each  target 
will  he  listed.  Since  there  are  usually  8  eigenrays  per  target,  and  128  targets 
per  target  depth,  it  follows  that  1024  lines  may  be  printed.  This  increases  the 
running  time  of  the  program  substantially. 

A  problem  often  encountered  in  sonar-system  analysis  is  to  determine  the 
variation  of  the  signal -to-noise  ratio  with  target  depth.  Since  reverberation  is 
independent  of  the  target  depth,  repetition  of  the  rather  lengthy  reverberation 
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computation  may  be  suppressed  by  punching  TOTAL  in  columns  41-50.  If, 
instead,  VOLUME  is  punched  in  columns  41-50,  surface  and  bottom  reverbera¬ 
tion  will  be  computed,  and  only  volume  reverberation  will  be  suppressed.  An 
example  of  a  valid  process  card  follows: 


Process 


print 


VOLUME 


OOOCO  000000308000000030000000000  0Q0G0  GG  0  0  0  0  0  0  0  0  0  Ok 

1  :  3  1  S  6  !  8  S  1011 1!  I1MIMSIHI 112621  22!)2W!’J!>»»%}l)n>!<UJ()MI»<54l424)4<<S<CO<l<](tSIS2S3S<Sa 

1  1  1  1  111  111  I  I  111  11  1  1 1  L-U-- -  1  1  111  I  1  111  1  1  UUO^ 


INITIAL  INCLINATION  ANGLES 

Initial  inclination  angle  cards  specify  those  rays  that  are  to  be  plotted.  *  As 
indicatedin  table  3C,  the  first  four  fields  on  the  card  generate  the  angles  and  the. 
last  four  limit  the  number  of  vertexes  and  boundary  intersections.  A  positive 
integer  punched  in  columns  37-40  will  cause  data  to  be  printed  as  well.  For 
example,  the  following  card  will  plot  rays  having  initial  inclination  angles  of  -5, 
-4,  -3,  -2,  -1,  -0,  and  +1  deg  and  print  the  coordinates  of  every  fifth  ray  segment: 


ioi  oooaooouo  ooooocsoo  oaocaacooooooooooisooooooooccoootsiottcoiioK  onionoo 
1 1 1 1 1 ! 1 1 1  i  i  1 1 1 1 1 1 1 1 1  1 1 1 1 1 1 1 1 1 1 1 ii 1 1 1 1 1 1 1 1 ii m 1 1 1 m mi i  mini;  nnnin  i 


22223222222222222222222222  222 


122222222  222222222222 


22222222222 


A  ray  will  be  terminated  when  it  has  undergone  2  surface  intersections,  1  bot¬ 
tom  intersection,  10  upper  vertexes,  or  10  lower  vertexes. 

\ 

END  CARD 

The  last  card  of  the  input  deck  should  be  the  end  card: 


000000 0000 000 3fl GOO OC3003000Q 0000 OG0OOG 
i  ?  «  s  5  :  i  j  uni;  nu  ijinniu;!  2i  ni)  2<  H  »?;)»  is  is  in: 

1  I  1  I  I  1  I  I  I  I  1  I  1  I  I  1  I  I  I  I  I  I  I  I  I  1  I  I  I  I  I  I  I  I  I  I  I  1 


*These  rays  are  not  used  in  the  acoustic  eigenray  computation. 
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constructs  a  ray  diagram  using  the  indicated  initial  inclination  angles.  If  the 
inclination  angles  are  missing,  the  ray  plot  will  be  suppressed.  The  @  EOF 
card,  which  is  necessary  whether  or  not  rays  are  plotted,  returns  control  to  the 
main  program,  where  propagation  loss,  reverberation,  signal-to-noi3e  ratio, 
and  probability-of-detecticn  curves  are  constructed.  Subsequent  cases  are 
similar  to  case  1,  except  that  the  data  which  do  not  change  from  one  case  to  the 
next  need  not  be  repeated.  The  end  card  terminates  the  execution. 


RAY  DIAGRAMS,  EXAMPLE  1 

As  an  example  of  a  NISSM  n  run,  let  us  consider  the  construction  of  ray 
diagrams'  (example  1),  which  required  24  s  of  Univac  1108  execution  time.  A 
typical  input  deck  is  listed  in  table  4.  As  stated  earlier,  axis  cards  are  unnec¬ 
essary  if  default  values  are  satisfactory.  The  circles  in  figure  25  denote 
input-velocity  data  before  curvature-of-the -earth  corrections  were  performed. 
However,  the  continuous-gradient  approximation  to  the  velocity  profile,  indi¬ 
cated  by  the  solid  line  in  figure  25,  does  include  curvature  corrections.  This 
difference  becomes  more  apparent  in  deep  ocean  profiles.  We  see  from  table  5 
that  only  tour  approximating  functions  (five  breakpoints)  were  required  to  fit 
the  data. 


The  ray  diagram  in  figure  26  corresponds  to  a  sonar  placed  on  the  sound- 
channel  axis  at  76  m.  Rays  were  traced  every  0. 2  deg,  from  -7. 0  to  +7. 0  deg, 
and  were  terminated  upon  intersecting  a  boundary  surface  or  reaching  the  maxi¬ 
mum  range  of  5  km.  In  this  particular  example,  the  maximum  number  of  vertexes 
was  not  attained.  A  partial  listing  of  printed  output  appears  in  table  6. 

The  ray  diagram  in  figure  27  is  similar  to  that  in  figure  26,  except  that  the 
sonar  was  placed  at  40  m,  which  is  36  m  above  the  sound-channel  axis. 
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Table  4.  Input  Deck  for  Constructing  Ray  Diagrams,  Example  1 


A  XQT  NISSM2 
COMMENT 

EPSTEIN  PROFILE 

A  EOF 

KANGb  axis 
uEPTm  axis 
VELOCITY  AXIS 
MAXIMUM  kANGE 

SONAh  depth 
VELOCITY  PROFILE 
M  M/S 

0.000  1502.19783 

4.000  1502.12858 

10.000  1501.99547 

14.000  1501.88285 

22.000  1501.58461 

30.000  1501.18430 

38.000  1500.59710 

46.000  1499.88046 

54.000  1499.0ol91 

62.000  1498,^0630 

70.000  1497.89307 

74.000  1497.55792 

76.000  1497,53779 

78.000  1497.58156 

86.000  1497.94980 

94.000  149B.b20S5 

102.000  1499.95718 

110.000  1501.12684 

118.000  1502.16418 

126.000  1502.99650 

134.000  1503.82023 

142.000  1504.06651 

150.000  1509.37602 

158.000  1504 ,5b627 

166.000  1509.72716 

174,000  1504,82071 

180.000  1504.8o935 

A  EOF 
PROCESS 

-7.0  47.0  0.2 

Cl  EOF 

SONAk  DEPTH 
PROCtSS 

-7,0  +7,0  0.2 

A  EOF 
END 
A  FIN 
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KM 

M 

5.0 
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Figure  27.  Ray  Diagram  Corresponding  to  a  Sonar 
Placed  Above  the  Sound- Channel  Axis,  Example  1 
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PROPAGATION  LOSS  VERSUS  RANGE,  EXAMPLE  2 


Propagation  loss  versus  range  is  computed  whenever  a  target-depth  card 
appears  in  the  input  deck.  For  example,  the  cards  listed  in  table  7  will  gener¬ 
ate  the  velocity  profile  in  figure  28,  the  ray  diagram  in  figure  29,  and  the 
propagation  losses  at  the  two  target  depths  in  figures  30  and  31. 


Table  8  converts  the  temperature-salinity-depth  data  to  velocities,  which 
are  then  plotted  in  figure  28.  See  table  9  for  the  corresponding  continuous- 
gradient  parameters.  By  making  use  of  the  minimum-range  card  and  inserting 
a  new  range  axis,  we  are  able  to  enlarge  a  particular  region  of  interest,  in  this 
case,  the  first  convergence  zone. 


Example  2  required  37  s  of  execution  time. 


If  the  word  PRINT  were  punched  in  the  process  card,  one  would  obtain  the 
eigenray  printout  of  table  10.  The  problem  would  then  require  41  s  to  be 
processed. 
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Table  7.  Input  Deck  for  Computing  Propagation 
Loss  versus  Range,  Example  2 


&  X01  NISSM2 
COMMLNT 

TYRkmENIA.I  StA-JUNE  1969 

t.  eof 

hANGt  AXIS  °*°° 

UEPTm  AX.S  °*° 

VELOCITY  AXIS  l*1-*0 

VELOCITY  PROP  ILL 

►T  DEG  F  /1000  J 

O.o  71.8  3&.2a 

Jo.O  71.3  38.24 

25.0  70.8  38. 23 

50.0  70.7  3b. 20 

bS.O  70,8  3«.lo 

100.0  62.3  38.11 

150.0  60.2  38.  Oti 

175.0  68.0  33.00 

200.0  57.6  38,00 

300.  u  56.9  36.04 

400.0  57.2  36.29 

500.0  57.4  38, 3b 

1500.0  57.2  3o.62 

£600.0  56.5  36.55 

4000,0  56.0  3').  4B 

6000.0  56.0  38.40 

9640.0  56.0  38.40 

c,  EOA 

OOTfoM  UtPTH  984 

bottom  Loss 

FREGuENCT  5 

sonah  depth  J 

MAXIMUM  MANGE 

sonah  tilt  anGlc 
sonah  response 
LEG  DO 

-90.0  18.0 

-16.0  18,0 

-15.0  14.6 

-14.0  16.1 

-13.0  17.4 

-12.0  15.8 

-11.0  13.2 

-10.0  H.O 

-  9.0  8,8 

-  8.0  7,0 

-  7.0  5,2 

-  6.0  3.6 

-  5.0  2,3 

-  4.0  1.0 

-  3.0  0,5 

-  2.0  0.2 

-  1.0  0.1 

0.0  0.0 

2.0  0.2 


55.0  11.9 
10000.0  5.0 
1.56  6.0 


9840. 0  FT 
3.0  MGS 
3.55  KCPS 
lb.O  FT 
55.0  RYD 
5.0  DFG 
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Table  7  (Cont'd).  Input  Deck  for  Computing  Propagation 
Loss  versus  Range,  Example  2 


3.0 

0.4 

5.0 

1.0 

7.0 

2.4 

8.0 

3.3 

9.0 

4.7 

10. 0 

5.8 

11.0 

6.8 

12.0 

7.9 

13.0 

9.0 

14.0 

10.4 

15.0 

12.2 

16.0 

13.9 

17.0 

15.8 

18.0 

*8.0 

90.0 

18.0 

a  EOP 
pKOCLSS 

0.0  20.0  1.0 

a  EOh 

hAUGL  AXIS 

PKCPaGATiON  lOSs  axis 
MNI.XUM  kANGE 
TARGtT  Ot-PTH 
PROCESS 
a  EOK 

lAHGtT  OlPTM 

PKOClSS 
a  EOI 
lMD 
a  Fli 


OEG 

2 

2 

35.0 

55.0 

4.0 

KYD 

40.0 

1*0.0 

5.0 

UP 

35.0 

KYO 

25.0 

n 

150.0 

FT 

4 
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Table  8.  Conversion  of  'Temperature-Salinity-Depth 
to  Velocities  for  Example  2 


DEPTH 

TEMPtRATURE 

SALINITY 

VELOCITY 

FT 

DtG  f 

/1000 

KM/S 

1 

0.0 

71.8 

38.25 

1.53126 

2 

16.0 

71.3 

36.24 

1.53061 

3 

25.0 

70.8 

38.23 

1.52993 

4 

bC.O 

70.7 

38.20 

1.52987 

5 

Ob.0 

70.6 

38.16 

1.52976 

6 

iOO.O 

62.5 

3B.11 

1.51734 

7 

150.0 

60.2 

36,06 

1.51366 

8 

175.0 

58.0 

38.00 

1.50992 

9 

200.0 

57.6 

38.00 

1.50934 

10 

300,0 

56.9 

38.09 

1.50870 

11 

400.0 

57.2 

38.29 

1.50999 

12 

bOO.O 

57.4 

38. 38 

1.51096 

13 

1500.0 

57.2 

38.62 

1.51592 

14 

2500.0 

56.5 

38.55 

1.51963 

15 

400o. o 

56.0 

38.48 

1.52623 

16 

oCOO.O 

56.0 

36.40 

1.54661 

17 

9t)40.0 

56.0 

38.40 

1.55616 
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SONAR  DEPTH  =  16  ft 
TARGET  DEPTH  =  150  ft 
RANDOM-PHASE  ADDITION  ■ 
COHERENT-PHASE  ADDITION 


45 

RANGE  (kyd) 


Figure  31.  Propagation  Loss  for  a  150-ft 
Target  Depth,  Example  2 
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REVERBERATION,  EXAMPLE  3 


The  input  deck  in  table  11  will  produce  the  velocity-profile  curve,  ray 
diagram,  and  reverberation  plot  in  figures  32,  33,  and  34,  respectively,  in 
3  min,  11  s.  Table  12  lists  the  corresponding  continuous-gradient-velocity 
parameters. 


Table  11.  Input  Deck  for  Constructing 
Reverberation  Plot,  Example  3 


A  XQT  NI5SM2 
COMMlNT 

REVERBERATION  -  RUN  6  -  SACLANTCEN 


A  EOF 

range  AXIS 
UEPTlt  AXIS 
VELOCITY  AXIS 
axis 

REVERBERATION  axis 
MAXIMUM  RANGE 
MAXIMUM  TIME 
SONAR  uEhth 
bOTTOM  LOSS 
HIND  SPElO 

freuuency 

PULSE  LEhGTH 
TRANSMITTER  LEVEL 
TRANSMITTER  RESPONSE 
RECEIVER  RESPONSl 
MOPI/ONTmL  BEAMaIDTH 

bOTToM  scattering  constant 

SURFACE  SCAHERIrG  CONSTANT 
VOLUME  SlATTERINg  STRENGTH 


M 

OB/YU 

0. 

-105. 

150. 

-105. 

270. 

-105. 

290. 

-97. 

AlO. 

-105. 

395. 

-97.5 

940. 

-9b. 

soo. 

-102. 

bOO. 

-95. 

730. 

-90. 

790. 

-84. 

b50. 

-73. 

920. 

-75. 

1050. 

-83.5 

1100. 

-85.7 

1315. 

-85.7 

2000. 

-105, 

2500. 

-105. 

A  EOF 

VELOCITY 

PROFILE 

M 

M/S 

0.0 

1515.0 

14.0 

1515,7 

27.0 

1514,0 

50.0 

1506.4 

120.0 

1507,0 

250.0 

1510.2 

940.0 

1513.0 

2500.0 

1547. s 

0.0  50.0  10.0  KM 

0.0  2.5  5.0  KM 

1.S0  1.55  5.0  KM/5 

0.0  50.0  10.0  RFC 

-100.0  0.0  5.0  OB 

50,0  KM 

50.0  SEC 

150.0  M 

3.0  MGS 

11.0  KNOTS 

3.5  KCPS 

0.58  SEC 


11  o.O 
SONAR 

sonar 

DB 

0.2 

RAD 

-27.0  Dh/YU 
♦  2.0  DB/YU 


A  EOF 
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Table  11  (Cont'd).  Input  Deck  for  Constructing 
Reverberation  Plot,  Example  3 


SONAk  RtaPONSf 

12.3  32.0  32.6 

53. A  33.  32.6 

34.3  34.7  35. 
24,9  29.5  29.4 

30.8  30.7  30. 

26.9  26,0  26,7 

20.3  19.  17.9 

20.  21.2  21.6 

0.8  5,6  4.5 

.3  .3  .4 

3.6  4.5  5,8 

22,  20,o  19 . 2 

19.4  20.5  21.6 

25,  24.0  24.6 

27.0  28.2  28. 4 

27.3  28.3  29.5 

31.4  31.2  31. 
31. 1  31.2  31.3 
31.3 

FROCtSS 

U.O  20.0 

3  EOF 

tUO 

2  ru 


32.9 

32. 7 

33.2 

29.4 
24. 
2u.5 
1  7.2 

21.2 

3.7 
«  5 
7. 

10.4 

22.8 

24.5 
20.4 
30  ,4 

30.7 
31.2 


33.  33.1 

32.7  32.7 
35.1  34.9 
29.6  29.8 
28.  27.9 

26.4  26, 

16.9  17, 
20,3  18,3 

2.9  2.2 

.7  ,8 

0.2  9.8 

17.8  17. 

23.9  25. 

24.5  24.7 

26.  27.6 

31.2  31.7 

30. 5  30.3 
31.1  30,9 


0.5 


COMPACT 

33.2  33.3 

32.9  33.2 
34.  33.1 

30.2  30.4 

27.7  27.6 

25.5  24.8 

17.2  17.7 

15.7  12,9 

1.4  .7 

1.2  1.6 

11.5  13,5 

lb. A  16.9 
25.2  25.3 
25.2  25.8 
27.  26.5 

31.7  31.8 

30.4  30.5 

30.8  30.9 


33.4  33.5 

33.5  33.9 

32.2  31. 

30.6  30.6 

27.  27. 

23.5  22. 

18.2  18.9 

10.4  8.4 

.4  .3 

2.2  2.8 

17.  22.4 

17.5  18.4 
26.3  25.1 

26.5  27,2 
26.3  26.6 

31.7  31.6 

30.7  30.9 
30.9  31.1 


DEO 


3  3  3 


80 
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PROBABILITY  OF  DETECTION,  EXAMPLE  4 


In  this  final  example,  we  consider  the  effect  of  varying  array  tilt  angles  on 
the  probability  of  detection.  Table  13  lists  the  input  deck  from  which  figures  35 
through  44  are  plotted;  table  14  contains  the  continuous-gradient  parameters. 
Example  4  required  5  min,  23  s  of  computer  time. 
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Table  13.  Input  Deck  for  Computing  Probability 
of  Detection,  Example  4 


*  XOT  NlSSM2 
COMMlNT 

COMPARE  WITH  FAST  N1SSM 

ft  EOF 

KANGl  axis 
UEPTh  axis 
VELOCITY  AXIS 
TIME  AXls 

PKOPaGAT iON  LOSS  AXIS 
KEVEHBEKmTION  AXIS 
SIGNAL  Tu  NOISE  RATIO  axis 
bOTToM  OcPTH 
SONAR  DEPTH 
TARGET  DtPTH 
FHEQuENCY 
ftIND  SPLEO 

aTTEnUATiON  COEFFICIENT 

bOTTuM  Scattering  constant 
SURFACE  SCATTERING  CONSTANT 


FT 

OB/YU 

0.0 

-  70.0 

2s. 0 

-  o9, 1 

75.0 

-  56. S 

12s.  0 

*  56, b 

17b. 0 

-  69.1 

225.0 

-  80.0 

300.0 

-  80.0 

400.0 

-  90.0 

1790.0 

-  90.0 

1810.0 

-120.0 

12000.0 

-120.0 

A  EOF 

UOTToM  LOSS 

VELOCITY 

profile 

FT 

FT/S 

0  0 

4870.0 

100. 0 

48/1.0 

2«0.0 

4843.0 

bOO.O 

4810.0 

2000.0 

4708,0 

.1000. 0 

4795,0 

12000.0 

4946.0 

ft  EOF 

TRANSMITTER  LEVLc 
HORIZONTAL  UEAMkIDTh 
blREcTIVITY  INuLa 
PULSE  LENGTH 
AVCRaOINu  TIME 

probability  of  f«lse  alahm 
target  SIRENGTh 

uANUaIOTh 
CLOSING  SPEED 
MAXIMUM  range 


0.0 

70.0 

7.0 

KYD 

0.0 

12.0 

6.0 

KFT 

4.75 

4.9b 

S.O 

KFT/S 

0.0 

70.0 

7.0 

SEC 

60.0 

120.0 

6.0 

DB 

-80.0 

♦  80.0 

5.0 

UP 

-80.0 

♦  80. 0 

5.0 

l)B 

12000.0 

FT 

20.0 

FT 

300.0 

FT 

3.5 

KCPS 

12.0 

KNOTS 

30.0 

DEG  F 

-27.0 

DB/YD 

♦2.6 

DB/YD 

4.0  MGS 


148.0 

OB 

10.0 

DEG 

26.0 

DB 

0 

15.0 

100.0 

10.0 


0.5  SEC 
0.5  SEC 
0001  XI 
DB 
HZ 

KNOTS 
60.0  KYD 
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Table  13  (Coat'd).  Input  Deck  for  Computing  Probability 
of  Detection,  Example  4 


MAXIMUM  IIME 

80.0 

SEC 

RECEIVER  RESPONSE 

SONAR 

TRANSMITIER  RESPONSE 

SONAR 

SONAk  response 

COMPACT 

02.8  32.0  32.8  32.9 

33. 

33.1 

33.2 

33.3 

33.4 

33.5 

33.4  33.  32.8  32.7 

3?. 7 

32.7 

32.9 

33.2 

33.5 

33.9 

J4,3  34.7  3S.  3a, 2 

35.1 

34,9 

34. 

33.1 

32.2 

31. 

29.9  29. a  24.4  29.4 

29.6 

29.8 

30.2 

30.4 

30.6 

30.8 

30.8  30.7  30.  24. 

28. 

27,9 

27.7 

27.6 

27. 

27. 

28.9  26.0  2b. 7  20.5 

26.4 

26. 

25.5 

24.8 

21.5 

22. 

20.3  19.  17.9  17,2 

lo,9 

17. 

17.2 

17.7 

18.2 

18.9 

20.  21.2  21.6  21.2 

20.3 

18.3 

15.7 

12.9 

10.4 

8.4 

0.8  5.8  4,5  3.7 

2.9 

2.2 

1.4 

.7 

.4 

.3 

•3  ,3  .4  ,5 

.7 

.8 

l.E 

1.6 

2.2 

2.8 

3.6  4.5  5.8  7. 

8.2 

9.8 

11.5 

13.5 

17. 

22.4 

22.  20. o  19.2  lo.4 

17.8 

17. 

16.8 

16.9 

17.5 

18.4 

*9,4  20. a  21.6  22.8 

23.9 

25. 

25.2 

25.3 

25.3 

25.1 

25.  24.0  24.6  24,5 

29.5 

24.7 

25.2 

25.8 

26.5 

27.2 

27.8  28.2  28.4  2o,4 

28. 

27.6 

27. 

26.5 

26.3 

26.  b 

27.3  28. a  29.5  30.4 

31.2 

31.7 

31.7 

31.8 

31.7 

31.6 

31.4  31.2  31.  30.7 

30.5 

30.3 

30.4 

30.5 

30.7 

30.9 

31.1  31.2  31.3  31,2 
31.3 

31.1 

30.9 

30.8 

30.9 

30.9 

31.1 

IHANsMITcER  tili  angle. 

20.0 

DEG 

RECEIVER  TILT  Ai.oLE 

20.0 

deg 

SONAR  TJuT  ANGLE 
PROCESS 

20.0 

DEG 

-20.0  0.0 

A  EOP 

+0.5 

OEG 

2 

THANSWIT1ER  TILI  ANGLE 

10.0 

DEG 

RECeIVER  TILT  Ai.oLE 

10.0 

OEb 

SONAR  TluT  angle 
PROCeSS 

10.0 

DEG 

iil  jii!  ■&,  il^'ii 


SONAR  DEI 
TARGET  DC 
RANDOM-1 
COHERENT 


j 

t: 


1 

i 
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REMARKS 

The  last  example  illustrates  that  the  NISSM  II  input  deck  is  relatively  com¬ 
plicated  when  compared  with  those  of  other  ray-tracing  programs.  Although  this 
will  somewhat  inconvenience  the  user,  he  will  find  compensations  in  the  numer¬ 
ous  input-output  features  of  the  computer  program. 

We  have  purposely  omitted  comparisons  with  actual  data  in  order  to  keep 
this  document  unclassified.  The  interested  reader  should  see  references  21, 
22,  and  23,  where  the  predictions  for  the  Continuous  Gradient  Ray-Tracing 
System  (CONGRATS)15  are  superimposed  on  the  Mediterranean  convergence- 
zone  and  reverberation  measurements.  NISSM  II  and  CONGRATS  differ  basically 
with  respect  to  program  efficiency,  but  the  predictions  for  both  models  are 
similar  and  the  program  results  agree  reasonably  well  with  the  measurements. 

It  is  felt  that  NISSM  II  is  a  validated  program  and  should  be  used  in  evaluating 
active  sonar  systems,  providing  the  limitations  of  the  model  are  understood. 
For  example,  several  irregularities  appearing  in  the  volume-reverberation 
plot  in  figure  34  are  a  result  of  caustic  curves  intersecting  a  strong  scattering 
layer.  The  treatment  of  caustics  in  this  situation  is  questionable;  therefore,  the 
data  represented  by  the  irregularities  may  be  erroneous  or  at  least  smaller. 
Hopefully,  further  analysis  will  result  in  improvements  to  NISSM  n. 
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C;-  3^ 
5_"i 


if  g^2)  =  0; 


if  g^  ^  *0;  and 


».M 

gf-l) 


Bo  = 


(A- 13) 


% 


if  gf_1)-0. 


gf  =  2g(n‘1) 


(A- 14) 


3 

if 

iS 


i 


•  *  *  Once  the  g®  have  been  set,  tests  are  performed  to  ensure  that 

a.  Depths  of  relative  minimum  Vq^  are  not  adjacent  to  those  of  relative 
maximum  v^\  and 

/jO  2 


[«t  *«?  «r 


.4 

1 


1 

I 


If  test  (a.)  is  not  satisfied,  then  an  additional  point  (of  inflection)  is 
inserted  at 


s 

Vr 

1 


(i) '  _ 


z  = 


+ z<1+1> 

0  0 


(A- 15) 


and 


ft 

■s'- 
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If  test  (b. )  is  not  satisfied  and 
reduced  according  to 


g, 


(i) 


<  g 


(i) 


g 


(i+1) 


,  then  g^l+1^  is 


g 


(i+1)  _  0.9999  [g(i)] 


’0 


(A-18) 


If  test  (b.)  is  not  satisfied  and  Jg^j  >  |g^  j>|g^+1^| ,  then  g^  is  reduced 
according  to  1  ‘ 


g 


2 

(i)  _  0.9999  [g(l)] 

e(i+1) 
g0 


(A-19) 


FinaUy,  if  |g<l)|  >  Jg(l)|  and  |gj*+1)|  >  |g^|,  an  additional  point  (of  inflection) 
is  inserted  according  to  equations  (A-15)  and  (A-16),  and  g(pf  is  given  by  the 
smaller  value  given  by  equations  (A-18)  and  (A-19) . 

The  next  step  is  to  compute  .  It  can  be  shownA1  that  for  c(z)  and 
dc(z)/dz  to  be  continuous  at 


,  ^  r  (i+1)  (i)l  (i) 

1  +  LZ0  -ZoJg2  = 


g(i)± 


g 


(i) 


g(i+1)  g(i) 
&0  50 


.(i+1) 


(A-20) 


if  Sq+1^  *  °» 


g 


Ji) 


1  + 


rz(i+i)-z{i)ie(i)=-^~ 

L  0  0  J  g2  (i) 


(A-21) 


if  gQ  =  0.  Only  positive  roots  are  acceptable,  and,  if  tests  (a.)  and  (b.l  are 
met,  at  least  one  root  will  be  positive.  As  the  last  step,  the  parameter  gjp  is 
given  by 


Al 

H.  Weinberg,  A  Curve  Fitting  Technique  for  Acoustic  Ray  Analyses, 
NUSL  Report  805,  14  March  1967. 
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g(i)  = 
S1 


gw  j 

i  r  ,h-i>  .  z(i)i 

r  ljo  \ . 

S2  j 

f  -  4" 

„(i+D  _  J 

l) 
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(A-22) 


In  practical  applications,  it  is  unnecessary  to  fit  the  input  data  exactly. 
Thus  the  number  of  breakpoints  in  the  interpolation  function  can  be  reduced  with 
a  corresponding  saving  in  computer  execution  time.  The  procedure,  which  con¬ 
tinues  until  j  equals  n,  is  outlined  below;  but  first  the  parameters  Vq  and 
gy,  i=l,2, . . .  ,n,  are  assumed  to  have  already  been  computed  and  modified 
that  tests  (a.)  and  (b.)  are  satisfied: 

•  The  index  i  is  initially  set  to  1. 


so 


•  The  program  determines  the  largest  index  j  and  parameters  and 
gjjp  such  that  a  function  m  the  form  of  equation  (A-l)  can  fit  the  data 
specified  at  zy  and  zy  while,  simultaneously,  the  maximum  devia¬ 
tion  between  the  input  velocities  and  the  values  of  the  fitting  function  at 
the  intermediate  depths  z^1+  \  zq+2K  •••»  does  not  exceed  0.2 m/s. 

•  The  intermediate  data  are  discarded  and  the  index  i  is  incremented  to  j. 
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Appendix  B 

EVALUATION  OF  RAY-TRACING  INTEGRALS 


RANGE  INTEGRAL 

Upon  substituting  the  velocity-interpolation  function  into  the  range  integral, 


Ax  = 


(B-l) 


one  obtains  an  expression  that  can  be  evaluated  using  elementary  rules  of  inte¬ 
gration.  However,  further  algebraic  manipulation  is  necessary  to  reduce  com¬ 
puter-truncation  errors.  The  results  are  summarized  below;  for  notational 
convenience,  the  superscript  (i)  has  been  dropped  and  z^  is  assumed  greater 
than  za. 


Define  the  constants  of  integration  a,  fi,  y  ,  p,  q,  and  r,  respectively,  by 

(B-2) 


a  =  v  -  C2  , 
o  v 


^  =  “g2+2g0' 


y  =  “go +  g 


i  * 


P=7-g2P  =  gl-2g0g2‘ 


q  ■=  -  /32  =  ar  -  g^  , 


(B-3) 

(B-4) 

(B-5) 

(B-6) 


r  -  gl  “  gQg2  * 


(B-7) 


and  the  variable  Y  by 


Y  =  o  +  2/3A  7.  +  tAz 


(B-8) 


B-l 


va 


3 
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Then  if 


-J 


dz 


(Y) 


1/2’ 


B1 

it  is  shown  by  Weinberg  that 


a 


f"  1  +  S2Az  2D  ! 

Ax=J  ~rv2dz=r  !1+s2(Az>av)+p 

•*  o  (Y)  v  * 


[y(2° "  il  ’ 


(B-9) 


(B-10) 


where 


<Az>  =  -~  (Az  +  Az  )  , 
N  'av  2  '  a  b 


(B-ll) 


D  = 


z,  -  z 
b  a 


1/2  1/2  ’ 
(Ya}  +  (Yb) 


(B-12) 


and 


if  |  yD2|  <  10-4  ; 


(B-13) 


H  =  2(-7f 1/2  tan-1  [(-7)1/2d]. 


(B-14) 


if  7D2  >  10  ^  and  7  <  0;  and 


H  =  2(7)  tanh 


1/2  ^-lp/2D] 


(B-15) 


if  J7D2  >10  4  and  7>0. 
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TIME  INTEGRAL 

Upon  substituting  the  velocity  function  into  the  time  integral, 


AT  = 


z 

/ 


c  |dz| 
v 1  1 


a 


4 


2  2 

z  c(z)  V  c„  “  c  (z) 


one  obtains 


AT  =  v  c  Ax  +  J , 
o  v 


where  Ax  is  given  by  equation  (B-10)  and 


J  = 


fb  Az^  +  g^z) 


J 


(l+g0Az)(Y) 


1/2 


dz . 


Further  analysis  shows  that  if  *  0,  then 


J  = 


cv  Ax  -  2pH  +  rl 


where 


J  a  + 


dz 


z&  (l  +  g2Az)(Y) 


1/2 


=  (r)~1/2  tanh"1 


1/2  <Yb>1/2  (i  +  »A°a)  -  (V1/2  (?  E0  +  ^  J 


1/2/1 
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(B-16) 


(B-17) 


(B-18) 


(B-19) 


(r) 


-r(Ya)1/2  (Yb)1/2  +(|  g0  +  pAaY|  g„  ♦  pAzJ 


(B-20) 


if  r  >  0,  and 
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r-M2 i 
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& 


•Jg 
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I  =  (-r)~1//2  tan" 


1/2  (2  B'o  +  pAza 

fs0+I,azb) 

r  -hVi/2  (vb)i/2 ♦< 

ho+paza) 

(5«0+pteb). 

(B-21) 


if  r  <  0.  There  is  no  need  to  evaluate  I  when  r  =  0  since  the  product  rl, 
which  appears  in  equation  (£-19),  approaches  zero  as  r  approaches  zero. 


RANGE  DERIVATIVE 


The  range  derivative  dx/dcv|z=Zb  is  obtained  by  differentiating  Ax,  as 
given  by  equation  (B-10) ,  with  respect  to  the  vertex  velocity  cv.  Since  p  is 
independent  of  cv, 


dx 

dc 


z=z. 


\D  dc  c  / 


2D  f  a<Az>av 

+r{h~K- 


+  p 


(B-22) 


We  now  consider  two  cases.  If  the  ray  does  not  vertex  in  the  layer,  za 
and  z^  axe  independent  of  cv. 


(B-23) 


90  _  da 

9c  g2  9c 
v  v 


lz=  g2  is. 

9c  °2  dc 

V  V 


+  g2Aza> 


2  da 

dc  ’ 
v 


(B-24) 

(B-25) 

(B-26) 


B-4 


JM**k*tf  *“ 


dc  ’ 
v 


d<Az> 


p_  p  AVacv  t  *V*cv\ 

"v  2  ((Ya)l/2  +  (Yb)l/2)  y  (Ya)l/2  (Yb)l/2 J  ’ 
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(B-27) 


(B-28) 


(B-29) 
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1  -yD2 

1  2  L 

v\J 

\2  dcy  Ddcv)  30^  L7\2D  7J 


(B-30) 


However,  if  the  ray  does  vertex  at  zb,  then  zb  will  be  a  function  of  Cy 
and  Yb  will  vanish  identically.  Equations  (B-22)  through  (B-26)  and  (B-30) 
are  still  valid,  but  equations  (B-27),  (B-28),  and  (B-29)  must  be  replaced  by 


d\ 

3o-=  0  ’ 

V 


d<Az>  dz 
av  1 _ b 

dc  2  dc  * 

V  V 


ap  1  D  dYa. 
*cv"“2Ya  dcv 


(B-31) 


(B-32) 


(B-33) 
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respectively,  where 


<X^2AV  «. 

2(yAzb+/3)  6cv 
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Appendix  C 

AMOS  SURFACE  DUCT  EQUATIONS0 


The  AMOS  equations'^  can  be  expressed  in  terms  of  the  following  variables: 

R,  horizontal  range  (kyd) 

Z^,  surface  layer  depth  (ft) 

Z  .  source  depth  (ft) 

Z  ,  target  depth  (1».) 

v 

f,  frequency  (kHz). 

When  the  sound-speed  structure  of  a  surface  duct  is  adequately  represented 
by  a  linear  sound-speed  gradient  y  and  an  average  sound  speed  C,  then,  in  any 
consistent  units,  the  range  between  successive  surface  reflections  of  the  limiting 
ray  is  approximately  2^2C/|-yi)ZL.  Hence  in  the  typical  situation  where 
7  =  0.018  s"1  and  C  =  4900  ft/s,  and  if  R  is  interpreted  as  the  number  of 
kiloyards  of  range  and  Zj,  is  interpreted  as  the  number  of  feet  of  layer  depth, 
r  is  twice  the  number  of  surface  reflections  of  the  limiting  ray  out  to  range  R, 
where 


xv 


Related  scaled  variables  introduced  for  convenience  are 

z  =  /Z  /Z~  ; 
x  ->/  x  L 


vJYV 


This  appendix  was  excerpted  from  references  cited  in  footnotes  C2  and  C3, 
C2 

H.  W.  Marsh  and  M.  Schulkin,  Report  on  the  Status  of  Project  AMOS 
Acoustic,  Meteorological,  and  Oceanographic  Survey) ,  1  January  1953  -  31 


December  1954,  NUSL  Report  255A,  23  March  1955. 

C3 

W,  H.  Wat3on  andR.  W.  McGirr,  An  Active  Sonar  Performance  Pre 
diction  Model.  NUC  Technical  Publication  286,  April  1972. 
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FIRST-ORDER  SURFACE  REFLECTION  AND  SHADOW  ZONE 

When  energy  has  been  reflected  at  least  once  off  the  surface,  corresponding 
to  <  r  <  r^  +  1/2,  the  following  equation  is  used: 

H  =  20  log(1000R/R  )+  «R  +  2(r-r  )  K  (z  ,z  ) 

v  1  t  X 

+  [I  -  2(r-ri>]  G  (zt,zx)  +  Nco  ,  (C-4) 


where 

K(zt, zx)  =  0.4Cf  [l0Zx  +  10Zt  +  10(Zt"Zx)]  ,  (C-5) 

with 


1,  f  <  8  (kHz) 

1/3 

(f/8)  ,  f  >  8  (kHz)  , 

and  Nco  is  the  low-frequency-cutoff  term  described  below.  For  this  zone  use 
either  equation  (C-3)  or  (C-4),  whichever  yields  the  lower  propagation  loss. 


SECOND-  (OR  HIGHER)  ORDER  SURFACE-REFLECTION  ZONE 

When  energy  has  been  reflected  at  least  twice  off  the  surface,  corresponding 
to  r  >r1  +  1/2,  the  following  equation  is  used: 

H  =  10  log(1000R/R  )  +«R  +  K(z  ,z  )  +  a  [R  -  (z~(r.  +  1/2)] 

v/  Z  X  S  V  Ju  X 

+  10  logflOOO^Z^  +  l/2)/R0]  +  Nco  ,  (C-6) 

where  ag  is  the  surface-reflection  loss  coefficient  in  decibels  per  kiloyard.  For 
this  zone  use  equation  (C-3)  or  (C-6),  whichever  yields  the  lower  propagation 
loss. 


C-3 
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SURFACE-REFLECTION  LOSS  COEFFICIENT 
The  expression  for  the  surface-reflection-loss  coefficient  ag  is  given  by 

<*g  =  T/Rb  (dB/kyd)  ,  (C-7) 

where  T  is  the  surface-reflection  loss  in  decibels  and  Rjj  is  the  bounce  dis¬ 
tance  corresponding  to  the  range  between  surface  contacts  for  the  surface-duct 
limiting  ray.  According  to  findings  of  Marsh  and  Schulkin,  C4  r  is  a  function 
of  the  product  of  frequency  and  wave  height,  as  given  by 


1.59J5T, 

fh  >  4.2691 

(C-8) 

10  log[l  + 

(fh/4.14)4],  fh<  4.2691, 

(C-9) 

where  f  is  the  frequency  in  kilohertz  and  h  is  the  mean  crest-to-trough  wave 
height  in  feet.  Figure  C-l  illustrates  surface  loss  per  bounce  as  a  function  of 
the  product  of  wave  height  and  frequency. 


1.0  10  100 
WAVE  HEIGHT  x  FREQUENCY  PARAMETER,  fh  (kHz  x  ft) 


Figure  C-l.  Surface- Reflection-Loss  Curve 


v*ffggL> 
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The  bounce  distance.  Rjj  is  determined  from 

Rb  =  JCL  ~  Cl/l500yQ  M  ’  (C“10) 

where  Cl  and  Cs  are  the  sound  velocities  (ft/s)  at  the  layer  depth  and  the  sur¬ 
face,  respectively,  and  yQ  js  the  sound-velocity  gradient  (s“l)  in  the  surface 
duct. 


LOW-FREQUENCY-CUTOFF  TERM 

Arase  and  others^S-Cll  have  reported  the  existence  of  low-frequency  cutoff 
in  surface-duct  propagation.  Thus  an  additive  cutoff-loss  term  has  been  incor¬ 
porated  that  is  based  on  an  approximation  to  the  normal-mode  surface -duct  model 


C5 

T.  Arase,  "Some  Characteristics  of  Long  Range  Explosive  Sound  Propa¬ 
gation,  "  Journal  of  the  Acoustical  Society  of  America,  vol.  31,  no.  5,  May  1969, 
pp.  588-595. 

C6 

M.  A.  Pedersen  and  D.  F.  Gordon,  "Normal-Mode  Theory  Applied  to 
Short-Range  Propagation  in  an  Underwater  Acoustic  Surface  Duct, "  Journal  of 
the  Acoustical  Society  of  America,  voi.  37,  no.  1,  January  1965,  pp.  105-118. 

C7 

M.  A.  Pedersen  and  D.  F.  Gordon,  Normal  Mode  Approach  to  Under¬ 
water  Sound  Propagation  (U),  NEL  Report  1407,  27  September  1967  (CONFI¬ 
DENTIAL). 

C8 

C.  S.  Clay,  "Sound  Transmission  in  a  Half  Channel  and  Surface  Duct, " 
Meteorology  International  Incorporated,  Project  M-153,  Technical  Note  Two 
2:1-17,  August  1968. 

C9 

W.  H.  Furry,  Theory  of  Characteristic  Functions  in  Problems  of 
Anomalous  Propagation,  Massachusetts  Institute  of  Technology  Laboratory 
Report  680,  1945;  and  Methods  of  Calculating  Characteristic  Values  for  Bilinear 
M  Curves,  Massachusetts  Institute  of  Technology  Report  795,  1946. 

CIO 

H.  W.  Marsh,  Theory  of  the  Anomalous  Propagation  of  Acoustic  Waves 
in  the  Ocean,  NUSL  Report  111,  12  May  1950. 

Cll 

D.  F.  Gordon  and  R.  F.  Hosner,  "Underwater  Sound  Propagation  ii 
the  Three-Layer  Ducts  Computed  by  Normal  Modes,"  Naval  Undersea  Warfare 
Center  Technical  Paper  107,  December  1968. 
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of  Pedersen  and  Gordon.  <^1^»  The  n-th- mode  amplitude  An  is  closely  approx¬ 

imated  by 


A„=U»(Zx^L»Un<VZL»e'100°TnR 


(0-11) 


where  Un(Zx/ZL)  and  Un(Zt/Zjj)  are  depth-dependent  functions  of  source  and 
target  depth  and  rn  is  a  mode-damping  factor  given  by 

v  K)1/3  **o*BwB  ■ 

where  Im(Mxn)  is  the  imaginary  part  of  the  n-th  eigenvalue.  Assuming  that  the 
empirical  AMOS  equations  correspond  to  the  fully  ducted  normal-mode  situation, 
that  Un(Zx/ZL)  Un(Zt/ZL)  is  a  slowly  varying  function  of  frequencv  near  cutoff, 
and  that  contributions  from  second-  (and  higher)  order  modes  are  negligible,  the 
cutoff  term  can  be  expressed  as 


N  =20  log(e'100°7’lR)  (dB)  . 


(C-12) 


Figure  C-2  shows  the  dependence  of  Im(Mx1)  on  M  and  p.  The  parameters  M 
and  p  are  defined  as 


M  =  <8'2tV/3V°s  ' 


(C-13) 


=  -  Vyi 


(C-14) 


The  quantity  is  the  below-layer  sound-velocity  gradient  (s_1).  Cutoff  effects 
start  near  M  =  2,3  and  become  extreme  near  M  =  1,0. C14  A  single  curve  cor¬ 
responding  to  p  =  -0.48  is  used  in  the  model.  Thus  in  computing  NCo  the  value 
for  Im(Mxi)  is  obtained  by  computing  M  and  then  interpolating  from  a  curve  of 
Im(Mx^)  versus  M. 


"Pedersen  and  Gordon,  op.  cit,. ,  footnote  C6. 
^Pedersen  and  Gordon,  op.  cit. ,  footnote  C7. 
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diction  System  (SHARPS)  program.  The  results  shown  in  figure  C-3  are  normal¬ 
ized  with  respect  to  the  FNWC  maximum  values  to  better  display  the  relative 
cutoff  effects.  This  comparison  shows  that  the  FNWC  cutoff  is  more  gradual 
than  the  NISSM  cutoff.  However,  this  difference  is  expected,  since  the  FNWC 
cutoff  uses  an  asymptotic  approximation  to  Im(Mxn)  that  underestimates  mode 
damping. 


